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Experiments on Porous-Wall Cooling and 


Flow Separation Control in a Supersonic Nozzle 


LEON GREEN, JR.,* anp KENNETH L. NALL** 


Summary 


Control of flow separation by fluid injection at one diverging 
boundary of a two-dimensional, transparent-walled de Laval 
nozzle was investigated by spark schlieren photography of dry 
nitrogen flows expanded from two stagnation temperatures 
(ambient and 500°F.) and several stagnation pressure levels. 
Injection conditions at the permeable boundary were varied by 
the use of three grades of porous stainless steel with nominal pore 
diameters of 10, 20, and 30 microns, through which nitrogen was 
forced by coolant reservoir pressures of 25, 50, and 100 psig, in 
addition to the case of no forced injection. Pressure distribu- 
tion measurements were made along the nonpermeable diverging 
boundary. It was found that flow separation at expansion ratios 
approaching the opiimum value for maximum thrust coefficient 
could be induced at the porous wall by a local injection mass 
velocity of the order of a few per cent of the local main-stream 
mass velocity. Separation at the solid boundary was not 
noticeably influenced by injection at the opposite wall, and the 
asymmetrical separation thus effected jet deflections of up to 10 
degrees at the lower stagnation-pressure levels. Variation of the 
wall heat-transfer condition by changing the stagnation tempera- 
ture did not significantly influence separation behavior. Tem- 
perature measurements at the reservoir face of the porous sec- 
tion, together with use of published correlations and of the Rube- 
sin analysis for estimation of stream-side Stanton Numbers under 
noninjection and injection conditions, respectively, permitted 
heat-transfer calculations which indicated that the effectiveness 
of the transpiration technique in controlling nozzle wall tempera- 
tures derives primarily from intimate fluid-solid contact in a 
porous material of high specific surface. 


Symbols 
A = area of channel cross section 
A* = throat area 
a@ = speed of sound 
b = channel width 
cy = friction coefficient 
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25, 1959. 
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Cp = specific heat of coolant gas 
Eu = Euler Number 

h = heat-transfer coefficient 

K = mixing-length constant 

k = thermal conductivity 

L_ = thickness of porous wall 

M = molecular weight 

M = Mach Number 

m = exponent 

m = mass flow rate of main stream 


m 
p = pressure 
Pr = Prandtl Number 

gas constant 

ry = temperature recovery factor 
Re, = length Reynolds Number 


s = distance from throat, measured along wall 
St = Stanton Number 
YT = temperature 
u = stream velocity 
v = normal injection velocity 
x = distance from throat, measured axially 
y = distance from manifold face of porous wall 
@ = viscous resistance coefficient of porous wall 
8 = inertial resistance coefficient of porous wall 
== viscosity 
p = density 

Subscripts 


= stagnation conditions 
= main-stream conditions 
= coolant 

= stream-side surface 
coolant-side surface 

= recovery 

= porous wall 

= zero injection 


Introduction 


ie THE DEVELOPMENT of long-duration rocket engines, 
the problem of nozzle cooling has been prominent. 
The available approaches for effecting controlled, 
steady-state wall temperatures may be categorized as 
convective cooling, in which heat transferred to the 
wall is removed by a coolant flowing through internal 


; 
ig & 4) 

4 

| 
a 

4 

> 

: 

4 


690 JOURNAL OF THE AERO/SPACE SCIENCES—NOVEMBER, 


passages; film cooling, in which coolant fluid streams 
injected at discrete locations coalesce to form a thin 
protective film adjacent to the wall; and porous-wall 
effusion cooling (also called ‘“‘sweat” or ‘‘transpiration”’ 
cooling) in which a permeable wall material is employed 
and an effectively continuous fluid injection distribu- 
tion is achieved. The cooling action of the effusion 
method, which is now generally agreed to be the most 
efficient in terms of degree of cooling effected by a given 
expenditure of coolant, derives from three factors: 
efficient transfer of heat from the wall to the infiltrating 
coolant is achieved because of intimate fluid-solid con- 
tact in a porous medium of high specific surface; if 
still in liquid phase on emerging from the wall, the 
fluid forms a cool insulating layer adjacent to the sur- 
face; the local ejection of mass into the boundary 
layer alters the velocity profile in a manner which 
decreases local heat (and momentum) transfer from 
the stream. 

It is of historical interest to note that the first re- 
ported experiments with the effusion or transpiration 
technique concerned nozzle cooling—in 1940 Meyer- 
Hartwig! successfully employed porous copper nozzles 
subjected to stream temperatures of 4,500°F. Porous- 
wall cooling of rocket combustion chambers for air- 
craft had been suggested by Goddard? in 1938, and its 
application to nozzle cooling was independently pro- 
posed by Skoglund’ in 1942. A brief series of ad hoc 
experiments was conducted by Jakob and Fieldhouse 
in 1943, using a small disc of porous brass subjected to a 
parallel flow of hot air, but the results were not reported 
until 1949.4 The first comprehensive investigation of 
porous-wall cooling was independently initiated in 
1944 by Duwez and his associates*: ® at the California 
Institute of Technology. In 1946, Moore and Grooten- 
huis’ proposed the technique for cooling gas-turbine 
blades. Since that time the effusion technique has 
become the subject of a voluminous and rapidly growing 
literature. It has been considered for application to 
flight structures (both for cooling of hypersonic vehicles 
and for heating to prevent ice formation on leading 
edges) and to the unit operations of chemical engineer- 
ing. The field has been reviewed by Grootenhuis and 
Moore’ for the period up to 1950, to 1952 by Eckert® 
and Broadwell and Sherman,’ to 1955 by Berger,'! 
through 1957 by Eckert, ef al.,!*.'* and into 1958 by 
the present authors.'‘ Experimental studies published 
as of this writing have employed bounded channel 
flows and unbounded flows along plates, cones, and 
blades. Analytical studies, which, with a few notable 
exceptions, generally deal with the more tractable case 
of laminar flow, have considered the above-mentioned 
situations, plus others, including constant-area channel 
flows with nonzero Euler Number, flows around blunt 
bodies, and flows with chemical reaction. No recent 
study of the nozzle problem appears to have been re- 
ported in the open literature. 

Although continued investigation of nozzle cooling 
by porous-wall effusion was thus desired, the principal 
objective of the present study was to determine the 
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feasibility of achieving controlled nozzle flow separa- 
tion. Boundary-layer control by distributed suriace 
suction had long been recognized as an effective means 
of delaying or preventing flow separation from an air- 
foil, and it was desired to evaluate the effectiveness of 
the inverse technique of surface injection in promoting 
jet detachment from the walls of an exhaust nozzle. 
Interest in the nozzle flow separation problem is due to 
the loss in nozzle thrust coefficient resulting from over- 
expansion of the flow to pressures below that of the 
surrounding atmosphere, an overexpansion which sub- 
jects single-stage rocket-propelled vehicles with nozzles 
designed for efficient high-altitude operation to a 
substantial performance penalty during take-off and 
low-altitude operation. 

The nozzle flow separation phenomenon has been the 
subject of several experimental and analytical investi- 
gations, which have been reviewed by the present 
authors.!® The experiments may be classified in two 
groups: those employing the hot combusion products 
from rocket motors and those using air or nitrogen 
expanded from ambient stagnation temperature to low 
static temperatures in the exit section. It was found 
that the results of nozzle pressure-distribution measure- 
ments from the various “hot flow’ experiments could 
be unified by an empirical correlation in terms of 
chamber pressure, atmospheric pressure, and pressure 
at the point of separation, independent of the gas prop- 
erties and nozzle divergence angle.'® '° The ‘‘cold- 
flow”’ studies, however, yielded results which conflicted 
concerning the divergence-angle effect and at best 
deviated slightly from the above-mentioned correlation. 
It was thought that this difference in behavior might 
be attributed to the different heat-transfer condition at 
the nozzle wall prevailing in the two types of tests— 
i.e., the stabilizing influence of a large heat flux from 
the stream to the wall might mask other effects. In 
view of the above situation, it was desired to study 
nozzle flow separation behavior under controlled stag- 
nation conditions, and, particularly, to investigate the 
effectiveness of boundary-layer injection in promoting 
incipient separation. 


Design of Experiment 


For this purpose, a two-dimensional, transparent- 
walled, de Laval nozzle similar to that of McKenney” 
was designed. Transparent walls of glass and later 
fused quartz were used. To promote smooth flow 
conditions at the throat, a converging nozzle entrance 
section was used to provide a continuous transition 
from the circular plenum cross section to the rectangular 
nozzle cross section. The entrance section was de- 


signed to be water-cooled if necessary, but this feature 
was not needed during the program described here. 
The nozzle flow-channel width was 0.5 in., the throat 
height, 0.2 in., the divergence half-angle 15°, and the 
exit expansion ratio 10. The nozzle-exit blocks were 
fabricated from stainless steel and were positioned by 
dowel pins. 


One block was tapped for static-pressure 
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Fic. 1. Schematic diagram of two-dimensional nozzle with 
porous-wall injection (porous section actually was located in up- 
per block). 


measurements along the diverging exit section. The 
other was provided with a section of porous stainless 
steel through which cold nitrogen could be injected 
from a manifold, as indicated schematically in Fig. 1. 
The assembled nozzle section is shown in Fig. 2. 

Previous investigations indicated that the measure- 
ment of the stream-side temperature of a porous wall 
was a source of difficulty, because the installation of 
thermocouple stations at this surface involved some 
perturbation of the permeability and coolant-flow 
distribution, and hence of the local heat-flow pattern. 
Because it was desired to use three easily exchangeable 
injection blocks employing porous sections of different 
permeability (each with five thermocouple stations) 
and because any perturbation of the coolant injection 
pattern was undesirable in that it might influence local 
separation behavior in a complicated manner, it was 
decided to measure the temperature of the manifold 
side of the porous-wall section and calculate the stream- 
side temperature. The injection nozzle blocks were, 
therefore, provided with ports through which spring- 
loaded, iron-constantan thermocouples could press 
firmly against the coolant-side porous surface. This 
nethod of porous-surface temperature measurement had 
previously been demonstrated by Grootenhuis, Mack- 
worth, and Saunders. '* 

Three porous-walled nozzle blocks were fabricated, 
employing (6.12-in. thick permeable sections of Micro- 
Metallic Type 316 porous stainless steel media with the 
properties tabulated in Table 1. The viscous and 
inertial-resistance coefficients, a and 8, which charac- 
terize the flow resistance of porous media over an ex- 
tended range of Reynolds Numbers according to the 
telation!® 


A(p?)/L = a(2RT/M)ppude + B(2RT/M) (pu)? (1) 
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TABLE 1 
Surface Viscous Inertial 
Nominal roughness resistance resistance 
pore dia. (RMS coefficient, coefficient, 
(microns ) microin. ) (in. ~?) (in. ~') 
10 160 7.41 & 105 2.25 105 
20 165 2.48 1.08 


30 210 0.42 0.51 


were determined from isothermal-flow calibrations of 
the individual sections, as shown in Fig. 3. 

Both the primary working fluid and the secondary 
injected fluid used in the present experiments were dry 
nitrogen gas from a liquid nitrogen supply. The flow 
rates and pressure levels were controlled by a system 
of Grove regulators and Fischer-Porter Flowrators 
capable of metering flows up to 1.0 lbm./sec. The 
nozzle-entrance stagnation temperature was controlled 
by varying the power input to the heat-transfer section, 
in which the gas flowed radially inward through a D.C. 
resistance-heated, porous graphite cylinder, as shown 
schematically in Fig. 4. The flow and temperature 
control systems, together with the manometer bank 
used to display the pressure distribution in the nozzle, 
the spark schlieren apparatus used to visualize the 
nozzle flow field, and the associated photographic and 
recording equipment, are described in detail else- 
where.!*. 


Experimental Conditions 


Systematic flow-separation experiments were per- 
formed at two main-stream stagnation temperatures, 
ambient (about 40°F.) and 5CO°F. <A few tests were 
completed using a stagnation temperature of 990°F., 
but fissure of one of the quartz walls (initiated during 
the 500°F. runs) under combined thermal and mechan- 
ical stresses prevented completion of a systematic 
series at this level. Nominal stagnation pressures of 
250, 300, 350, and 400 psig were employed in the 500°F. 
tests, but flow system capacity limited the pressure 
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Fic. 3. Flow resistance characteristics of porous wall media. 


levels in the ambient-temperature tests to the range 
200, 250, and 300 psig. Four injection conditions were 
employed, with coolant at ambient temperature forced 
by manifold pressures of 25, 50, and 100 psig in addition 
to the case of no forced injection (manifold sealed). 
Coolant-side porous-wall temperatures were measured 
only in the tests with heated main-stream flows. 

The pressure and velocity conditions prevailing in 
the nozzle are indicated by the profiles presented in 
Fig. 5, which were computed from the known area-ratio 
variation by standard one-dimensional gas-dynamic 
methods assuming isentropic flow and neglecting the 
small mass addition through the porous wall section. 
The Euler Number of the flow was calculated from 
the pressure and velocity profiles of Fig. 6 and is shown 
in Fig. 7. The apparent discrepancy between the 
results from differentiation of the pressure profile and 
those from the slope m = Eu of the velocity profile on 
logarithmic coordinates is attributed to the fact that 
m is the average slope of a flat but slightly inflected 
curve. The length Reynolds Number distribution for 
the various stagnation conditions employed is shown 
in Fig. 7. 


Experimental Results 


Schlieren Observations 


The nozzle flow fields observed by spark schlieren 
photography for flows expanded from ambient stagna- 
tion temperature are typified by the examples shown in 
Fig. 8 for the particular case of the medium-grained 
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(20 micron nominal pore dia.) wall section and 4 
nominal stagnation pressure of 300 psig. Fig. 9 pre. 
sents similar photographs of flows expanded from 
stagnation conditions of 500°F. and 350 psig. In each 
figure, the nozzle locations corresponding to the throat, 
the point of optimum expansion to one atmosphere 
pressure, the point of jet detachment from the porous 
wall, and the detachment point on the nonporous wall 
predicted by the empirical correlation of separation- 
pressure data'® are indicated. The separation point 
predictions appear to be more accurate in the case of 
the heated flows than for the unheated flows. It may 
be noted that even in the case of no forced injection 
some asymmetry of the separation pattern is evident. 
This behavior is attributable to the strong streamwise 
pressure gradient in the flow. Hence, gas forced out- 
ward through the porous wall in a region near the 
nozzle throat re-emerges downstream in a region of 
lower static pressure, thus producing some injection 
effect. 

The results shown in Figs. 8 and 9 indicate that as the 
injection velocity is increased, the point of jet detach- 
ment from the porous wall approaches the optimum 
separation point—i.e., that where the local static 
pressure equals | atm. Within the range of conditions 
studied, the effectiveness of injection was not noticeably 
influenced by the stagnation temperature. It should 
also be noted that the rate of fluid injection required to 
promote early separation is not great; the injection 
ratios at the point of jet detachment varied over a 
range from about 0.5 per cent to about 3 per cent in the 
88 different combinations of stagnation conditions, 
injection pressure, and wall material employed in the 
investigation.'!: As indicated in Fig. 10, integration 
of the injection velocity distribution—calculated from 
the nozzle pressure distribution by Eq. (1)—1from the 
leading edge of porous insert (s = 0.38 in.) up to the 
point of separation (s,.,) revealed that the total coolant 
expenditure amounted to only a few per cent of the 
main-stream flow rate, and no effort had been made to 
minimize this fraction. 


TABLE 2 
Experimental Conditions Under Which Typical Nozzle Flow 
Fields Shown in Figs. 8 and 9 Were Obtained 
(Porous wall media: 29-micron nominal pore diameter; 
atmospheric pressure 14.4 psia) 


Fig. 8 Fig. 9 


stagnation temperature (°F.) 48 500 
stagnation pressure (psig) 300 350 
main stream flow rate (lbm./sec. ) 0.68 0.59 
coolant temperature (°F.) 43 44 
coolant manifold pressure (psig) 
A * * 
B 25 25 
c 50 59 
D 100 100 
local injection ratio (pudw/pit:) at point 
of separation from porous wall (per 
cent) 
A * * 
B 0.51 0.68 
Cc 0.83 1.01 
D 1.21 1.5 


* Coolant manifold sealed. 
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Fic. 4. Schematic diagram of heat-transfer section. 
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Fic. 5. Distribution of gas-dynamic parameters in nozzle (neg- 
lecting mass addition). 
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Fic. 6. Distribution of pressure gradient, velocity gradient, and 
Euler Number in expansion section. 
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Fic. 7. Reynolds Number distribution in nozzle exit section for 
various stagnation conditions. 
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Fic. 8. Schlieren photographs of typical nozzle flows expanded 
from ambient temperature. 


55 ADDITION 


CALCULATED POINT 
PRL by 


POINT (NEGLECT 
MASS ADDIT 10M) 


ATED POINT 
PRESSURE 


ATED POINT SHORE 
om anes 


SE PARAT 
POINT (NEGLECT 
ADDITION) 


T rs 1 


Fic. 9. Schlieren photographs of typical nozzle flows expanded 
from 500°F. 
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Nozzle expansion ratios at separation vs. total injection 
ratio. 


Fic. 10. 


From the above semiquantitative studies, it is con- 
cluded that the porous-wall injection technique holds 
promise as a method of achieving controlled flow sep- 
aration in supersonic exhaust nozzles. Also noteworthy 
is the angle of jet deflection effected by the asymmetrical 
separation pattern. Deflections up to 10 degrees were 
observed in tests at the lower stagnation pressure levels. 
It appears possible that the use of locally-variable 
boundary injection from segmented manifolds with 
pressures regulated by a valve system might constitute 
a means of thrust vector control under altitude condi- 
tions permitting flow overexpansion. Such an arrange- 
ment would enjoy advantages over the conventional 
methods of control by chamber or nozzle gimballing or 
by the use of vanes or jetevators*! in that the actuating 
forces and inertia of moving parts would be low, and no 
moving parts would be subjected to heating and ero- 
sion by the stream. 


Static Pressure Measurements 


The pressure distribution along the nonpermeable 
diverging nozzle boundary corresponding to the photo- 
graphs of Fig. 9 are shown in Fig. 11, which also indi- 
cates the location of the optically-observed points of 
natural separation from the nonporous walls (as well 
as those of forced separation from the porous wall) for 
comparison with the estimated separation point. 
The results, which are typical of those observed under 
other conditions, indicate that the estimated separation 
station agress well with the points of minimum pres- 
sure, whereas the flow actually detaches from the wall 
slightly further downstream in a region of adverse 
pressure gradient. 


Heat-Transfer Calculations 


It was mentioned earlier that measurements of the 
temperature of the coolant manifold side of the porous 
wall were substituted for more difficult measurements 
of the stream-side temperature. The situation involved 
in calculation of the stream-side temperature was shown 
schematically in Fig. 1. The quantities which are known 
or measured directly are the porous-wall permeability 
characteristics; the nozzle geometry, A* and A(x)/A*; 
the stream mass flow rate, '; the stagna ion pressure 
and temperature, Py) and 7»; the coolant pressure and 
temperature, ~, and 7; and individual values of the 
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manifold-side porous-wall temperature, 7),,, at five 
stations. The quantities which may then be directly 
calculated are the distributions of mass velocity, Mach 
Number, static pressure, and static temperature in 
the main stream, and, hence, length Reynolds Number 
Re, based on fluid properties evaluated at main-stream 
conditions. In addition, the distribution of coolant 
injection mass velocity through the wall may be cal- 
culated from a knowledge of the distribution of pres. 
sure-square difference p,.2 — ),° across it because the 
temperature difference across the wall is small, and 
(in so far as fluid properties entering the injection-rate 
calculation is concerned) no significant error is intro- 
duced by considering it to be isothermal at known 
temperature 77,,. 

However, the stream-side wall temperature 7, 
cannot be precisely calculated without knowledge of the 
temperature recovery factor, r, and the heat-transfer 
coefficient, h. Thus, there is required a solution of the 
heat-conduction problem in the wall simultaneous with 
an evaluation of the effect of injection upon r and h. 
In view of the relatively fine texture and high specific 
surface of even the coarsest porous-wall medium, it 
may be safely assumed that the difference in tempera- 
ture between the porous solid and the gas permeating 
it is negligible compared to the total temperature levels 
involved. The heat conduction in the wall (assumed 
to take place only in the high-conductivity solid phase) 
thus is governed by the boundary-value problem 


ky (d?T,/dy*) — pu®uCp(dT,/dy) = 0 (2) 
T(0) = Tx, (3a) 
ky —| = A(T, — Ty,) (3b) 
dy 


Using average values for the thermal properties k,, and 
cy in the wall, the solution to this problem is 


(Tw; — Tw) /(T; — Ty,) = [1 — 
(4) 
Of the various analyses of transpiration cooling which 


have been reported, that of Rubesin®? was selected for 
evaluation of the effect of injection upon h and r, owing 
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to its success in correlation of experimental data on heat 
and momentum transfer in turbulent, compressible 
boundary layers in supersonic as well as subsonic flows 
with injection and with positive Euler Numbers. In 
view of the large Reynolds Numbers involved (10° to 
10’), the roughness of the porous walls, the destabilizing 
effect of injection, and the results of conventional nozzle 
heat-transfer experiments by Saunders and Calder™ and 
Pefley,** the conclusion that the boundary layer on the 
porous section would be turbulent appeared safe. 
Although Rubesin’s method involves many steps, it was 
adaptable to an iterative calculation procedure which 
could be programed for machine solution on an Elec- 
trodata E-102 digital computer. This procedure is 
described in detail elsewhere!‘ and will only be outlined 
here. 

The starting point of the calculation was to assume 
an initial value of 7\,, (7; = Tw. was used at first, but 
experience soon permitted a better first approximation 
and, hence, more rapid convergence of the iteration) in 
order to evaluate the various parameters involved in 
Rubesin’s analysis. With the distributions of stream 
mass velocity, injection mass velocity, and length 
Reynolds Number determined as described previously, 
the subsonic-flow data of Mickley and Davis” were ex- 
trapolated to higher injection ratios to obtain a starting 
approximation for the friction coefficient c,/2. This 
estimate provided the basis for an iterative calculation 
which converged upon a value of c,/2 satisfying all rela- 
tions for the initial estimate of 7», Although the basic 
analysis assumed that the Prandtl-Karman mixing- 
length constant A was fixed at a value of 0.392, un- 
affected by transpiration, the present calculation 
arbitrarily employed variable values of A again by 
extrapolating the results of Mickley and Davis to higher 
injection ratios. The validity of the process has not 
been assessed, but it appears as defensible as does the 
assumption of constant value of A in the region adja- 
cent to an effusion-cooled wall with large injection. 
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Fic. 13. Stanton Number vs. alternative injection parameter. 


The calculation of c,/2, however, constituted but one 
step in a larger iteration loop for determination of 77,,. 
Values of c,/2 from the internal iteration loop were used 
to evaluate Rubesin’s expressions for the Stanton Num- 
ber, taking the turbulent Prandtl Number as unity and 
evaluating the molecular Prandtl Number at local 
main-stream conditions. At this point, the external 
iteration loop was closed by calculation of 7),, from Eq. 
(4), and the process (including the internal iteration for 
c,/2) was repeated until convergence was achieved. 

The results of the Stanton Number calculations from 
all the heated-flow tests are presented in Fig. 12 in 
terms of the ratio St/Sts vs. Only 
data corresponding to stations within the nonseparated 
zone (as established by comparison with the schlieren 
photographs) are included. It was found that the 
injection ratios employed were so large that a long, non- 
linear extrapolation to the zero-injection condition 
would be required to establish values of St, from the 
experimental data. Therefore, Stx was estimated 
from the expression 


Stx = 0.0285 Re,~'/°Pr-?4 (5) 


where Re, is the Reynolds Number based on the dis- 
tance downstream from the throat measured along the 
wall. This type of approximation has been found 
satisfactory for correlation of conventional nozzle heat- 
transfer data.** *4 The calculations indicate that the 
injection ratios employed were sufficient to reduce the 
heat-transfer ratios by more than one order of magni- 
tude. It is seen that, despite the drastic mixing-length 
extrapolation involved in the modified application of 
Rubesin’s method, the St/St, points appear to fall 
within the band calculated by Rubesin® to cover the 
range of variablesO < 8,1 < 7,/7; < 3,and 10°< 
Re < 10°. An alternative presentation of the Stanton 
Number calculations in terms of an alternative injec- 


tion parameter, Re,, is shown in Fig. 13. 
This parameter is of significance in analyses of laminar- 
flow situations where it is proportional to the Blasius 
function evaluated at the wall. It is used here only for 
correlative purposes because it permits a linear extra- 
polation of the data in the range of low-to-moderate 
injection ratios which is of greatest practical interest. 
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parameter. 


The wall-temperature ratios (7, — 7T.)/(To — Tv) 
obtained from the calculations are presented as a func- 
tion of (p,%»/pit1)/Stx in Fig. 14. It is seen that the 
effect of injection is to reduce the stream-side surface 
temperature by more than an order of magnitude com- 
pared to that which would be experienced under adia- 
batic recovery conditions. The apparent divergence 
of the points is magnified by the logarithmic ordinate 
scale (which was necessary to distinguish the small 
temperature-ratio values which otherwise clustered close 
to zero) and is not regarded as significant. 

Corresponding values of the friction coefficient and 
temperature recovery factor derived from the calcula- 
tion procedure are shown in Figs. 15 and 16. The re- 
covery factor is seen to increase slightly with increasing 
injection ratios. This trend is surprising in view of the 
results of other (laminar-flow) analyses of single-com- 
ponent boundary-layer injection. However, Rubesin 
advises that his analysis probably does not yield pre- 
cisely valid results for the influence of injection on re- 
covery factor, and that the most positive statement 
possible is that the effect is small, with an approximate 
value of r = 0.9 being representative for the regime of 
low-to-moderate injection rates which is of practical 
interest. 


Concluding Recommendations 


The experiments described in the foregoing discussion 
indicate that the technique of effusion cooling by injec- 
tion of a gaseous coolant through a porous wall is highly 
efficient under conditions of strong flow expansion (high 
positive Euler Number). In addition, it was demon- 
strated that, when the injection technique is applied to 
one diverging wall of a two-dimensional, supersonic de 
Laval nozzle, early separation of the flow from the 
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porous wall may be promoted under expansion condi- 
tions approaching the optimum. By means of such 
controlled, asymmetrical flow separation, controlled jet 
deflections up to 10° (depending upon the chamber 
pressure and injection ratio) could be effected. Since 
the experiments were not designed to investigate this 
behavior over a wide range of conditions, it is recom- 
mended that the attractive possibility of using porous- 
wall injection techniques for control of the magnitude 
and direction of the nozzle thrust vector be evaluated 
experimentally by small-scale rocket thrust chamber 
tests using axially-symmetric nozzles with segmented 
injection manifolds. 


Effect of Film-Coefficient and Recovery Factor 
Assumptions 

As a means of evaluating the sensitivity of the heat- 
transfer calculations summarized above to the arbi- 
trary assumptions required, a second set of solutions to 
Eq. (4) was performed subject to the assumption that 
the film coefficient and recovery factor are not affected 
by injection—i.e., St = St» and r = Pr'’3, The wall- 
temperature ratios so derived are presented in Fig. 17 
for comparison with the correlating line from Fig. 14 
and with the line estimated by Rubesin” for the case of 
conventional cooling of 100 per cent effectiveness—.e., 
when the coolant achieves full wall temperature. It is 
seen that the data reduced in this fashion fall somewhat 
below Rubesin’s estimate in the range of low injection 
rates. 

Fig. 17 indicates that the effectiveness of the effusion 
cooling technique derives principally from the efficient 
heat transfer between the porous solid and the infiltrat- 
ing fluid during its passage through the wall, a con- 
clusion supporting that of Lifton.” The additional 
beneficial effect of mass addition to the stream upon the 
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boundary-layer conditions governing the film coefficient 
is thus illustrated by the first-order difference between 
the two lower curves of Fig. 17. It would thus appear 
that differences in prediction of the effect of mass 
transfer on the film coefficient by the various analyses 
currently available may be considered to be of second- 
order importance in the design of porous-wall cooling 
systems. In cases where the transfer coefficient per se 
(rather than the wall-temperature reduction effected) is 
of primary concern, as in the case of skin-friction esti- 
mates, these differences between the various analytical 
predictions, of course, become significant. 
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An Experimental Study of Jet-Flap 
Compressor Blades' 


EDWARD L. CLARK, JR.,* ann DONALD EARL ORDWAY** 


Cornell University 


Summary 


The results of a preliminary experimental investigation to 
determine the feasibility of using the jet flap to improve the sec- 
tion characteristics of an axial-flow compressor blade are pre- 
sented and discussed. The airfoil tested was the NACA 65- 
(12)10 blower blade with a thickened trailing edge. Internal 
design of the blade is described and details of the resulting jet 
flow are given. Also included are wind-tunnel design and test 
procedures for the two-dimensional cascade used in the test. 

Test results are presented in the form of the measured turning 
angle, pressure rise, and lift coefficient. They are examined with 
particular reference to the prevention of rotating stall. 


Symbols 


= blade span 
blade chord 


Cr = force coefficient 

C; = jet-momentum coefficient 

C; = section lift coefficient 

F = force 

g = coordinate parallel to tangential direction 

Zo = tangential spacing between blades 

M; = jet momentum 

p = static pressure 

Pa = atmospheric pressure 

Pr = total pressure 

pe; = total pressure in the jet 

q = dynamic pressure 

V = flow velocity 

y = coordinate perpendicular to downstream flow direc- 
tion 

vo = defined as gp cos Bo 

Ap/q: = nondimensional static-pressure-rise parameter 

a = angle between upstream flow direction and blade 
chord 

B = angle between flow direction and cascade axis 

6 = angle between jet direction and blade chord 

9 = flow turning angle, (8: — B2) 

p = mass density of air 

o = solidity; ratio of blade chord to blade tangential 
spacing 

Subscripts 

a@ = component in axial direction 

m = referred to vector-mean velocity 

¢ = component in tangential direction 

1 = upstream of blade row 

2 = downstream of blade row 
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Superscript 


* = referred to equivalent downstream velocity 


Introduction 


QO of the most serious problems involved in the 
design of high performance axial-flow compressors 
is the problem of rotating stall. The phenomenon of 
rotating stall has been widely discussed in the liter- 
ature! * and is given only briefly here. The stalling 
occurs when a row of blades, such as the rotor or stator 
of a compressor, has a high angle of attack relative to 
the upstream flow, as shown in Fig. 1. This frequently 
happens during rapid acceleration of a turbojet engine— 
for example, prior to take-off or during in-flight throttle 
adjustments. A small disturbance can then increase 
the flow angle of one of the blades to the point where 
it will stall. This stalled condition may create a re- 
gion of separated flow which will deflect the upstream 
flow. The next blade will then have a higher effective 
angle of attack and it may also stall. At the same 
time, the deflected flow decreases the angle of attack 
of the blade below and tends to prevent or remove the 
stalled condition. In this manner, the stall region 
“rotates” around the blade row. 

One approach to the elimination of rotating stall 
would be the use of boundary-layer control to delay 
flow separation and, hence, increase the stall angle of 


the blade. However, while boundary-layer control 
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Fic. 2. Cascade test section as viewed from tunnel exit. 


has been extensively studied and tested on isolated 
airfoils, only a few attempts have been made to apply 
it to the blades of turbomachines. Three of these 
which have been examined are: (1) suction through 
slots,? (2) suction through porous walls,‘ and (3) 
slotted blades. None of these references give data 
which can be related to the specific problem of rotating 
stall, although Schultz’s experiment‘ showed that there 
was an increase in lift for an angle of attack less than 
the stall angle. 

Another possibility was suggested to the authors 
by Prof. F. E. Marble of the California Institute of 
Technology. His idea was to utilize the jet flap which 
has been successfully tested on isolated airfoils, for 
example, see reference 6. This consists of a high-speed 
stream of air blown from the trailing edge of an airfoil, 
usually at an angle to the chord. The jet of air, 
similar to a solid flap, increases the airfoil circulation. 
Furthermore, it reduces separation effects by adding 
high energy air to a region which would normally have 
a low velocity. Therefore, these effects should allow 
a jet-flap blade in cascade to operate at a lower angle of 
attack than a conventional blade with the same load. 
They should also give a higher stall angle if the jet 
deflection angle is not too large.’ 


Description of the Cascade Tunnel 


It was of prime interest to isolate airfoil character- 
istics, such as turning angle and lift coefficient. There- 
fore, the airfoils were tested in a two-dimensional cas- 
cade tunnel rather than in an actual compressor. 

The Cornell 10 X 25 in. low-speed cascade tunnel is 
of open circuit design.’ It utilizes a large settling 
chamber with eight fine-mesh screens to provide uniform 
flow at the test section. No method was employed for 
removing the test section wall boundary layer although 
this would have been desirable. 

The inlet air angle 8; was 45° and the solidity o was 
1.0. The cascade was composed of seven airfoils. 
Each had a 4 in. chord and an aspect ratio of 2!/.. 
Although it would have been advantageous to have had 
all airfoils of the jet flap type, tunnel construction 
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necessitated a conventional airfoil at either end of the 
cascade. The test section end walls were cambered 
and could be adjusted to match the cascade. 

Compressed air for the jet flaps was supplied from 
two manifolds to both ends of the five “blowing’’ air- 
foils. An overall view of the cascade test section is 
shown in Fig. 2. 


Description of Airfoils 


The airfoils selected for this test were NACA 65- 
(12)10 compressor blades with trailing edges thickened 
to a 2 per cent radius.* This airfoil section was used 
since it has a trailing edge suitable for the internal con- 
struction desired and is representative of compressor 
blades having well-known characteristics. 

To avoid three-dimensional effects, it was necessary 
that the internal duct provide a constant exit jet ve- 
locity along the span. In an effort to achieve this 
result, the internal channel converged from a relatively 
large plenum chamber to a small parallel slot. This 
spanwise slot was 1/2 in. long and 0.013 in. thick. As 
shown in Fig. 3, the slot was inclined at an angle of 
15.9° to the chord. Preliminary tests showed that it 
was essential for the exit slot to be of uniform thickness 
across the span. Therefore, the airfoils were cast of 
Meehanite to avoid residual stresses associated with 
materials such as cold-rolled steel. For structural 
strength, 1/10 in. thick internal ribs were placed at 
one-third and two-thirds of the span. 

The spanwise total-pressure distribution, as measured 
by a probe inserted one-tenth of an inch into the slot, 
was fairly uniform as shown in (a) of Fig. 4. From 
this measurement, it was assumed that the jet effects 
except at the rib stations were essentially two-dimen- 
sional along the span. However, it is to be noted that 
the total-pressure distribution one chord downstream 
of the trailing edge was not as uniform, see (b) of Fig. 4 
for the test airfoil. An earlier “test block’’ used in 
developing the internal configuration showed better 
results, (b) in Fig. 4. It also had internal ribs but a 
lower ratio of rib to slot thickness than the test airfoil. 
Limited experiments indicate that the downstream 
uniformity of the jet sheet depends critically on this 
ratio, among other variables. 


Test Procedures and Measurements 


Initial tests at the design angle of attack indicated 
that the airfoil characteristics became practically inde- 
pendent of Reynolds Number at a value of 270,000, 


Fic. 3. Cross section of the NACA 65-(12)10 compressor 
blade, with 2 per cent trailing-edge radius, showing internal 
duct design. 
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Fic. 4. Spanwise total-pressure (in in. of Hg) distribution with 
zero tunnel velocity and an airfoil plenum pressure of 20.4 in. of 
Hg: (a) inside exit slot and (b) one chord downstream of exit slot 
of airfoil (solid line) and ‘‘test block’’ (dotted line). 


which is typical for cascade tests." '' With a blade 
chord of 4 in., this allowed an upstream velocity of 140 
ft. per sec. that was low enough to avoid compressibility 
effects. The corresponding dynamic pressure was 
22.1 Ibs. per ft.”, and all tests were made at this pressure. 

The cambered end walls were adjusted for each air- 
foil angle-of-attack setting. This was necessary in 
order to satisfy as well as possible the criteria of in- 
finite blade row conditions. These criteria are: (1) 
the upstream and downstream flow directions should 
be uniform along the blade row, usually with the up- 
stream direction set parallel to the tunnel axis at one 
chord, and (2) the upstream static pressure should be 
uniform along the blade row. It was found that the 
end walls could be positioned to satisfy either of the 
conditions, but not both. An initial test showed, how- 
ever, that best results were obtained with uniform 
static pressure so this criterion was used. This re- 
quired, of course, that the upstream flow direction 
be measured and included in the turning angle. Be- 
cause of the two conventional end blades in the cas- 
cade, it was found that neither of the infinite blade 
row conditions could be rigorously satisfied with jet 
coefficients other than zero. The test procedure used 
was to set the end walls for optimum conditions with the 
jets off and then leave them in that position when the 
jets were on. 

The upstream flow direction was measured with a 
cobra-head yaw probe.'* The probe head was placed 
at mid-span of the airfoils, one chord upstream of the 
cascade. Readings were taken with it located on the 
tunnel centerline and then 2'/. in. on either side of this 
position. The upstream static pressure was measured 
with wall orifices located along the bottom cascade wall, 
one-half chord ahead of the blade row. At this posi- 
tion, any “‘local’”’ influences ot the airfoils were found 
to be negligible. Total pressure ahead of the cascade 
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was measured by two probes located at mid-height of 
the tunnel, | in. from the end walls. 

Downstream of the cascade, the flow direction was 
also measured with a cobra-head yaw probe. The 
probe was set at mid-span of the airfoils, one chord 
downstream of the cascade. All measurements were 
made in the airfoil passages adjoining the center airfoil, 
However, since the yaw probe gives erroneous results 
in the presence of a total pressure gradient, all measure- 
ments were taken outside of the airfoil wake regions. 
The total pressure distribution one chord downstream 
of the cascade was determined with a seven-tube non- 
integrating rake. 

As defined, the jet coefficient is the nondimensional 
momentum of the jet. Thus, the coefficient can be 
determined by direct measurement of the jet momen- 
tum. With zero tunnel velocity, the jet total pressure, 
and, hence, momentum, was measured with the total- 
head rake at mid-span of the airfoil. Measurements 
were made one-half chord and one chord downstream 
of the airfoil trailing edge. These results, practically 
equal, were averaged. The resultant momentum was 
converted to a nondimensional form by using the 
product of blade area and upstream dynamic pressure. 
Jet coefficient values used in the present test were ().15, 
().33, and 0.52. 


Data Reduction 


Calculation of the forces on a cascade was based on 
pressure and momentum changes as determined from 
measurements taken upstream and downstream of the 
cascade. As anticipated, the results were not in exact 
agreement with the results expected from an ideal, two- 
dimensional, infinite cascade. 

To begin with, the flow in the cascade test section 
was not actually two-dimensional. This was _prin- 
cipally due to the nonuniform jet sheet, secondary flow 
in the blade channels, and the interaction of wall and 
airfoil boundary layers at the airfoil tips. This last 
effect introduced a large error in measurements of the 
static-pressure rise through the cascade. Fig. 5, show- 
ing the theoretical pressure rise, gives an indication of 
the magnitude of this error. The measured value is 
less than the theoretical because the walls experience 
the same pressure rise as the airfoils and, in addition, 
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Variation of cascade-pressure-rise coefficient with turning 
angle for C; = 0. 
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Momentum control surface for a cascade of jet-flap 
airfoils. 


Fic. 6. 


have a well-developed boundary layer prior to the 
cascade. Hence, a thick boundary layer is produced, 
creating a smaller effective cross section downstream 
of the cascade. This ‘‘contraction’’ causes the exit 
velocity to be larger than that predicted by continuity 
relations based on a constant tunnel height. There- 
fore, a lower pressure is obtained. 

On the other hand, the jets issuing from the trailing 
edges of the airfoils entrain part of the main flow down- 
stream of the cascade. This lowers the velocity of the 
main flow there, and so increases the static pressure 
rise. In fact, subsequent examination of the test data 
will show that the measured pressure increment is 
greater than would be expected from the increase in the 
measured turning angle. Hence, it would be inter- 
esting to calculate the entrainment but this is difficult 
because of the flow acceleration by the exit ‘‘contrac- 
tion.” 

Finally, as it was previously mentioned, the tests were 
actually run with an inlet air angle which was as much 
as 4° greater (one chord upstream) than the geometric 
angle of 45°, since the upstream flow could not be set 
parallel to the tunnel axis. However, the turning 
angle and lift coefficient change only very slowly with 
inlet angle in this range. 

The force on a cascade is frequently resolved into 
axial and tangential components. To calculate these 
forces, consider a control surface containing one airfoil 
as shown in Fig. 6. The axial force F,, taken positive 
in the upstream direction, is given by the momentum 
equation, or 


£0 
F, = ef (V2 cos B2)*bdg — p( Vi cos B1)*bgo + 
0 
— pi)bgo (1) 
where the density p is assumed constant. Also, letting 
the tangential force F; be positive if the velocity change 


is in the usual direction, shown in Fig. 7, we have, 
similarly, 


F, = p Vi2bgo cos sin By — 
pb f V2? cos sin 62 dg (2) 


It is now convenient to define an “equivalent”? down- 
stream velocity V2* which will give the same resultant 
force but will be independent of tangential position. 
Assuming the flow direction in the wake to be the same 
as the flow direction outside of the wake and taking the 


equivalent velocity in this direction, we find 
(Vet)*= Vitae (3) 
0 


Now, far downstream of the cascade the static pressure 
is atmospheric,t so for incompressible flow, Eq. (3) 
may be rewritten, 


(V2*)? = (2/pgo cos = pa)dy (4) 


where the y axis is taken perpendicular to the flow 
direction—that is, yo = go cos #. The integrand 
(pi, — Pa) was determined by a survey with the total- 
head rake, and the integral evaluated graphically. 

The force coefficients Cry, and Cr, based on the 
upstream dynamic pressure g; and the blade area bc 
are then, 


Cr, = Fi/qbe = (2/c)[(V2*/ Vi)? cos? Be — 
cos” B; + (Ap/2q:)| (5) 


Cr, = F./qbc = (2/c) [cos sin — 
(V2*/ V1)? cos sin (6) 


where Ap = (p2 — pi) and o =c/qo, the cascade solidity. 
A standard procedure in cascade tests is to define the 
lift and drag as components of the resultant force per- 
pendicular and parallel, respectively, to the vector 
mean velocity V,,. Since this velocity is the vector 
average of V; and V2*, its position is determined by 


Bm = arctan [(Vi, + V2,*)/(Vo, + Va*)) (7) 


as shown in Fig. 7. Hence, the axial and tangential 
force coefficients may be resolved into lift and drag 
coefficients. 

The jet coefficient C; is defined as 


C; = (8) 


where the jet momentum ./; is given by 


M; = f pV°bdy (9) 
jet 
t In order to observe the pressure recovery due to the ‘‘ejector 
effect,’’ the sidewalls should be extended further downstream as 
practiced for highly loaded turbine blades. 


Fic. 7. Typical velocity diagram for an axial compressor. 
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Fic. 8. Variation of cascade static pressure rise with angle of 


attack for 8; = 45°, « = 1.0, and NACA 65-(12)10 blade section 
with 6 = 15.9°. 
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Fic. 9. Variation of turning angle with angle of attack for 
8; = 45°, o = 1.0, and NACA 65-(12)10 blade section with 
6 = 15.9°. 


Far downstream of the airfoil, the flow is essentially 
incompressible; therefore, 


C= f “(by = paddy (10) 
je 


The integral of Eq. (10) and, hence, Eq. (10) is evalu- 
ated from a total-pressure traverse of the jet with zero 
tunnel velocity. 


Results and Discussion 


The blade section characteristics are given in Figs. 
8-10. The present tunnel configuration did not afford 
sufficient accuracy to determine the drag coefficient. 

It can be seen that the stall angle was increased by 
the use of the jet flap. Since measurements near the stall 
condition were difficult to make with the present tunnel, 
the actual increase was not determined. However, the 
important result obtained was that at the maximum 
jet coefficient of 0.52 the lift was increased more than 
30 per cent throughout the angle-of-attack range. Or, 
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a comparable lift coefficient can be obtained at a much 
lower angle of attack. For example, if the blades were 
normally operated at a design angle of attack of 14°, 
they could now be used at an angle of 7° with a jet 
coefficient of 0.52 to obtain the same lift. Thus, 
through the increase in stall angle and increase in lift, 
a row of jet flap blades would have a greatly improved 
chance of avoiding rotating stall. 

The efficiency of the jet flap blade would be an im- 
portant factor if the blowing were used continuously, 
The most likely source of air for the jets would be a 
bleed from the later stages of the compressor. There- 
fore, with continuous blowing it is possible that the 
increased pressure rise due to the jets would be offset 
by the losses incurred in supplying the jet air. How- 
ever, if the jets were used only at critical periods, the 
benefits gained by prevention of rotating stall would 
outweigh any losses in efficiency. It should be noted, 
though, that the blade design does not take advantage 
of the increase in circulation at high speeds when the 
jet flap is not used at design speed. Therefore, all of 
the gain in angle-of-attack range discussed above 
cannot be realized. 

The present tests constituted only a preliminary in- 
vestigation and, hence, were crude in many ways. 
Future tests will include a porous wall cascade test 
section with flexible end walls and a single-stage, axial- 
flow jet-flap compressor. These tests, which will in- 
clude measurements of the loss coefficient and effi- 
ciency, will help to determine if the jet-flap blade can 
be a solution to the rotating-stall problem. 
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Fic. 10. Variation of lift coefficient with angle of attack for 
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6 = 15.9°. 
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The Equation of Energy Balance for 
Fluttering Systems With Some Applications 
in the Supersonic Regime 
JOHN D. C. CRISP* 


Massachusetts Institute of Technology 


Introduction 


[' IS WELL KNOWN that the development of a critical 
flutter oscillation can be attributed to the existence, 
at an appropriate air speed, of a right (characteristic) 
set of amplitude ratios and phase angles (eigenmodes) 
between certain component motions—being then neu- 
trally stable, the system is potentially capable of energy 
gain (or loss) as a result of some further perturbation. 
It seems natural then that a study of the eigenmodes, 
rather than the eigenvalues (critical flutter speed and 
frequency parameter) should lead to a clearer under- 
standing of the flutter mechanism and of those system 
parameters which are most significant in allowing it. 
This is indeed the case; and it will be shown how such 
a study leads to some useful tools of assessment for the 
practicing aeroelastician. The emphasis in analysis 
here passes from interest in the solution, by whatever 
means, of the characteristic flutter determinant for the 
eigenvalues to an assessment of the eigenmodes which 
is best made by consideration of the equation of energy 
balance. This, in turn, illuminates the equations of 
motion themselves so that a cursory inspection of the 
matrices of the flutter coefficients can be very reveal- 
ing as to the system behavior and stability. Such 
assessment as is possible with the concepts outlined 
below provides a useful guide in formulating and inter- 
preting conventional flutter analyses. 


These considerations are particularly pertinent today 
when the inevitable complexity of multi-degree of 
freedom systems tends to obscure the basic mecha- 
nisms, and the flutter engineer is in danger of losing 
his physical understanding of the problem. 


For a harmonically oscillating conservative system, 
the Integral of Energy is a first integral of the equations 
of motion and, in fact, is a statement of the law of con- 
servation of energy. A fluttering system, on the other 
hand, is subject, inter alia, to dissipative forces so that 
the integral of energy leads to an equation of balance 
relating the total energy of the system to the energy of 
the impressed aerodynamic forces or generalized forces, 
as they are customarily termed. Whereas the con- 


Presented at the Aeroelasticity Session, IAS 27th Annual 
Meeting, N.Y., Jan. 26-29, 1959. 

* Formerly, Research Engineer. Now, Research Assistant 
Professor of Applied Mechanics, Polytechnic Institute of 
Brooklyn. 


ventional analysis of the equations of motion leads to 
a complete statement of system behavior, the equation 
of energy balance, being a first integral of the equa- 
tions of motion, is less informative in this sense. It 
states, in essence, those amplitude and phase relations 
necessary for the existence of flutter when the system 
displacement is represented by a truncated modal 
series. It does, however, reveal those system proper- 
ties which are basic and significant in flutter in a gen- 
eral way which elucidates the results of the more pedes- 
trian parametric solutions of the equations of motion; 
and it is capable of development as a criterion of sta- 
bility which finds useful practical application. 


The equation of energy balance has been used in 
various restricted forms since Routh! applied it to sta- 
bility considerations of mechanical systems, a topic 
further elaborated by Rayleigh.* Their work formed 
the basis of its application, by Frazer and Duncan,’ 
to stability determinations for some rather special 
forms of two-dimensional binary flutter. Garrick‘ 
used a two-dimensional form of it to investigate pro- 
pulsive effects of a flutter-like wing motion, while 
Barton,*® Greidanus,’ Rott,® and Lambourne'! employed 
it to derive limiting reduced frequencies for flutter to 
exist in the case of two-dimensional wing bending- 
torsion flutter. Of these authors, Rott carried the ideas 
furthest, to a point, in fact, where he deduced a simple 
criterion of flutter stability. This was expressed in 
geometric terms as a “limit circle’ in the complex 
amplitude ratio plane which is a special case of the 
“energy circle’ derived independently by the present 
author'® from a criterion of “‘ultimate stability’’ for 
binary fluttering systems of any type of coupling, and 
two- or three-dimensional in nature. 


Two experimental applications of the energy balance 
equation are to be found in the ‘‘flutter engine” of Ellis 
and Duncan,’ and the wind-tunnel measurements of 
Halfman.!” 


The generalization of the statement of the equation 
of energy balance to any finite number of degrees of 
freedom of a three-dimensional system was elegantly 
made by Duncan,’ and he emphasized the importance 
of the ‘‘skew-symmetric aerodynamic stiffness com- 
ponents” of the generalized forces, which will be 
elaborated below. The development given here is 
quite general; it is applicable to any finite number of 
degrees of freedom of a three- (or two-) dimensional 
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linear system, caters for time-dependence of system 
parameters, is valid for divergent and subsident os- 
cillations. Specialization is then made to neutrally 
stable oscillation, thence to binary three-dimensional 
flutter. A criterion of ‘“‘ultimate stability’’ is estab- 
lished for this case. Several examples of supersonic 
flutter include wing-bending-torsion flutter, wing- 
aileron flutter, and panel flutter, and are designed to 
illustrate the application and limitations of the equa- 
tion of energy balance, both as a means of understand- 
ing the essential behavior of fluttering systems and as 
a design tool. 


Equations of Motion 


A familiar form of the equations of motion can be 
written in matrix form as 


= 0 (1) 


where the matrix of g is a column of » generalized co- 
ordinates and the elements 6, of the m X nm square ma- 
trix are 


= —[1 — (1 + ig) + Qi (2) 


Here a, are generalized masses, g; are modal damping 
coefficients, VW; = w,/w) the modal frequency ratio, 
while 2 = w/w and k = wh/U are the required eigen- 
values. The generalized (aerodynamic) forces Q;, are, 
of course, complex. It will be convenient at first to 
ignore structural damping (place g; = (0) and to rewrite 
Eqs. (1) and (2) in a general matrix form more amen- 
able to energy derivations, thus: 


Ag + Bg + Cq = 0 (3) 


where A, B, C are square matrices containing as ele- 
ments the generalized masses, the imaginary, and 
the real part of the generalized forces added to the 
structural stiffness forces, respectively. Thus, Q, = 
Cy. + 2ikb;, when the motion is sinusoidal. The g 
and its derivatives are to be understood, therefore, as 
column matrices. It will add little complication to 
generalize the situation to allow a time dependence of 
all system parameters so that Eq. (3) becomes 


A(t)g + + + Cg = 0 (4) 


It is worth while now relating the forces of Eq. (4) to 
the basic work functions from which they may in part 
be derived, thereby gaining a formal description of the 
system. 


The Force Potentials 


Quite formally, the forces implied by the matrices 
A, B, C are derivable in part from scalar functions 
or potentials 7, F, V, respectively. We define 


T = 
F = [q’B(t)q)/2 
V = [q’C(t)q)/2 (5) 


where 7 is the kinetic energy quadratic in the general- 
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ized velocities, F is Rayleigh’s dissipation function, and 
V is a “generalized potential energy.” The prime g’ 
etc., indicates the transpose of g, etc. Fora fluttering 
system, the matrix A is symmetric, and for the usual 
choice of generalized coordinates, diagonal. The ma- 
trix B, in the absence of any mechanical or structural 
damping, is of aerodynamic origin and usually non- 
symmetric; in particular cases it too may be symmetric, 
as, for example, when piston theory oscillatory forces 
are applicable, the only restriction being the need for 
profile symmetry about the mean profile chord. The 
matrix C is almost without exception nonsymmetric, 
and consists of aerodynamic and also structural 
components, the latter of which by themselves usually 
form a diagonal matrix. It is of vital interest to take 
note of the form of these matrices of elements because 
it reveals in a general way the importance of couplings 
as well as the probable behavior of the system. 

It is a great convenience to write the matrices B, C 
as a sum 


B= B, + Bi; C=G4+C (6) 


where B,, Ci, are symmetric and By, Cy are skew-sym- 
metric—i.e., obj, = or obj, = (bj — dyj)/2. It 
will transpire that the skew-symmetric components* 
of the ‘‘stiffness’’ matrix C are of great significance in 
flutter. 

The definitions (6) substituted in Eq. (5) lead, by a 
property of operations on skew-matrices, to a form of 
the potentials 7, F, V which involve only the symmetric 
A, B,, Ci; since x’Boxx = x'Cx = 0. The scalar V 
thus has as its coefficients both structural and (sym- 
metric) aerodynamic parts and is a “generalized” 
potential energy in this sense. 

It remains to recall that the coefficients A, B, C, A 
of Eq. (4) are in general functions of the reduced fre- 
quency k and Mach Number VM. This is not always 
so, however; for piston theory they are all independent 
of k. 


The Total Energy Gradient 
The total energy of the system is given'® by 


Ln = E (7) 
j= 


where /; are the generalized momenta 0L/0g,, and L is 
the generalized kinetic potential 7-V. Because of the 
nonconservative nature of a fluttering system (through 
the dissipation function F and the explicit time de- 
pendence of L), the total energy E is not a constant. 
Using an expression for g from Eq. (4) substituted in 
the relation 


dL /dt = ((L/24))4) + + 


it is straightforward to show that 


* Rayleigh? calls these gyrostatic terms. 
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THE EQUATION OF 


(d/dt)(E) = + — 
E=T+V 


for the special forms (5) of the potentials, where there 
arises a new generalized dissipation function 


6 = (1/2)g’(B + A/2)g = + A/2)g (10) 


which degenerates to F when the generalized masses 
are time independent—.e., A = 0. Thus, the rate of 
change with time of the total energy of the system is 
determined by the nature of the generalized dissipation, 
the magnitude of the skew-symmetric components of the 
aerodynamic ‘‘stiffness’’ forces, and the time variation of 
the symmetric components of those same forces. It is 
clear that the stability of the system is likewise de- 
pendent on the same factors, inasmuch as increasing/ 
decreasing E implies gain/loss of energy to the system. 
It is of interest to note that the time rate of change of 
the damping forces, B, does not appear in Eq. (9). In 
fact, for a fluttering system experiencing acceleration 
or undergoing structural stiffness changes due to kinetic 
heating it is the last term which primarily governs the 
transient behavior. For conventional flutter, C,; = 0, 
itis the second term involving C, which dominates and 
primarily determines system stability. This is because 
fluttering systems are characterized by large skew-sym- 
metric aerodynamic stiffness couplings and the system 
can often be assessed by a scrutiny of these. 

It is only in the single degree of freedom system that 
the dissipation function ® or F is the prime behavior 
determinant, since then the right-hand side of Eq. (9) 
becomes — + ay 2). For ay, = 
(y = 0 this implies the well-known result that stability 
is determined by the sign of the damping coefficient by). 

In the conventional flutter case, 4 = C, = 0, 


(d/dt)(E) = (d/dt)(T + V) = -—(2F + q’Cq) (11) 


which is equivalent to Duncan's result.* 


Energy and Stability 


Eqs. (9) and (11) indicate that system stability is 
assured if the expressions within the respective brackets 
of the right-hand side are always positive for the whole 
range of a system parameter (such as airspeed) of in- 
terest. (7° is necessarily positive definite, and we 
assume for the moment that so also is V’.) Looking 
first at Eq. (11), we note that for aerodynamic (flutter) 
systems the dissipation function is not necessarily 
positive definite (as is usually the case for mechanical 
systems) but is so in most cases. The quantity g’C.g, on 
the other hand, is as likely to be negative as positive, 
and it is the behavior of this quantity which almost 
invariably dominates system stability, and may, for 
example, vitiate all attempts to induce stability by 
artificial increase of damping—i.e., by algebraically 
increasing the value of the dissipation function. 

We may state these ideas more explicitly, on the 
basis that both 7 and V are positive definite, by writing 
G = thus: 
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(a) If G = 0, for which a sufficient condition is C, = 
0, i.e., if there are no skew-symmetric stiffness couplings, 
then, provided F is positive definite (or |B,|; > 0), 
the system is stable irrespective of the skew-symmetric 
damping couplings By. 

(b) If B, = 0, ie., if there are no skew-symmetric 
damping couplings, then, provided F is positive definite 
(or |B, 5 > 0), the system may be stable or unstable. 

Similar considerations apply to the more general 
Eq. (9) from which three conclusions can immediately 
be drawn. Firstly, we note the complete absence of the 
time-variation of the dissipative forces, B. Secondly, 
the time rate of change of the stiffness forces clearly 
dominates the transient behavior, and it is of interest 
that of these the symmetric components (C,) are the 
significant ones and, moreover, they act in a sense 
opposite to the skew-symmetric components Cy. 
Finally, time rate of change of generalized masses 
provides an effective additional damping effect, which 
may be positive or negative. An example is later given 
incorporating some of these phenomena. Suffice it to 
say here that they are readily interpretable physically. 
The variation + | A |, for example, may be due to time- 
variation of the density ratio yw induced either by 
altitude or mass changes. The variation + 'C,| could 
arise simply through accelerated flight or, additionally, 
through structural stiffness changes due to kinetic 
heating. For example, in the case that 4 = 0), tran- 
sient effects on stability are governed by the sign of 
q'Ciq. Positive discriminants of this expression assure 
its positive definiteness. For a binary system, ” = 2, 
for example, this requires 


— > 0; (én) > (12) 


If satisfied, Eq. (12) implies a destabilizing influence. 

It should be noted that the definition of instability 
as the increase of system total energy is an “‘instan- 
taneous’’ one in the sense that it ignores the system 
condition with indefinite increase of time, and, in fact, 
avoids the difficulty of interpreting stability in phys- 
ical terms. This leads to no ambiguity for linear 
systems with time-independent coefficients, but when 
these conditions are not met, the absorption of energy, 
if time-limited, may not necessarily result in unaccept- 
able system motion or catastrophic failure. For such 
systems, our definition of instability then strictly ap- 
plies only for the time interval of energy absorption, 
it still being necessary to interpret the resultant state 
of the system after the cease of energy-absorption as one 
of structural integrity or its lack. With this proviso 
the definition stands, and is used in the sequel. 

For aerodynamic systems, the generalized potential 
energy | need not be of positive definite form. In this 
case, the system may be subject to instabilities of a 
nonoscillatory nature, so-called divergence. Since 
flutter instabilities are of prime importance only in the 
absence of such divergence, we assume throughout that 
these do not exist so that the earlier assumption of 
positive definiteness of V is retained. 
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Fic. 1. Paraboloid surface of stability for binary flutter. 


The Equation of Energy Balance 


The integral with respect to time of Eq. (9) between 
limits f) and ¢, yields the gain or loss of total energy of 
the system over that time interval. Thus: 


aE = + Vis = 
ty 

-f [26 + — (13) 
to 


which is the equation of energy balance, relating the 
system total energy change to the external work. 
When the motion is periodic, and harmonic of circular 
frequency k, it is permissible to place 4 = 0, 4h = 
a 2r/k. It follows by virtue of the form of the poten- 
et tials T and V that for no transient effects (A, C inde- 


pendent of time), Eq. (13) reduces to 
AE = - f (2F + g’C.g)dt = 0 (14) 
0 


which is a necessary condition for flutter to exist. This 
expression depends only on certain of the aerodynamic 
forces, the modal amplitude ratios and phase angles, 
and involves no inertia or structural parameters ex- 
plicitly. 


Kinematically Admissible Flutter 
Configurations 


The implications of Eq. (14) become clearer upon the 
specification of the type of motion. For neutrally 
stable flutter oscillations, we may define the motion, in 
matrix form, as* 


[q;] = [Q, sin (kt + 7 = 1,2,... 


which upon substitution in Eq. (14) and after integra- 
tion yields 


,n (15) 


n n 


0 = AE = DY cos (6; — &) — 


Cy sin (e; | (16) 


* t here can be regarded as the nondimensional time Ut/b. 
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where the suffixes (1) and (2) refer to symmetric and 
skew-symmetric components, respectively, and is the 
reduced frequency. Define now the amplitude ratio 
R; and phase difference 


Ry = Q)/Qr; An = & (17) 


and the mean energy per cycle becomes 


e = = — > Rj Re cos Ax — 
1 


j=l j= 

sin Xx] = O (18) 
Eq. (18) is simply a statement of a condition that the 
n — 1 amplitude ratios R and the n — 1 phase angles 
A must satisfy for flutter to exist at a given reduced 
frequency k. It thus defines all kinematically admis- 
sible flutter configurations, which are determined by 
the symmetric (aerodynamic) damping forces and the 
skew-symmetric stiffness forces. The latter will, in 
general, be solely of aerodynamic origin, inasmuch as 
the matrix of structural stiffness force coefficients is 
symmetric (usually diagonal, in fact). In order to 
see how Eq. (18) can be developed toward practical 
ends, it is convenient to consider the simple binary 
flutter case. 


Binary Flutter and the Energy Circle 


In the event that only two degrees of freedom are 
coupled, n = 2, Eq. (18) is reducible to a simple form 
by placing Q, = Qe, « = 0 which normalization im- 
plies no loss of generality. Then R, = 1 and for sim- 
plicity let R; = R, Aw = «. Eq. (18) becomes 


= —bu{ + 2(x cos — Esin )R+ ¢} = 0 (19) 


where x = (bi + bn) 
= (C12 = /2byk 
= (20) 


This is the polar equation of a circle, such that a vector 
of modulus R and argument e¢ based on the origin of 
cartesian coordinates traces the circle of center (—y, £) 
and radius 7 given by 


(21) 


In other words, this circle defines kinematicaliy ad- 
missible flutter configurations of R, «. (There exist 
actually two further restrictions on admissible R, « 
pairs. First, because R is by definition positive, Eq. 
(19) indicates that Z = x cose — Esine <0. Secondly, 
in addition Eq. (19) considered as quadratic in R can 
have real positive roots only if Z2 — ¢ 2 0. Both in- 
equalities place bounds on the phase angle e but they 
are seldom of significant interest in relation to the 
energy circle.) Such are transitional between configu- 
rations lying within the circle (when e2 > 0) which per- 
mit energy gain/instability, and those without the 
circle (when e;2 < 0) which require energy loss/stability. 
Consider now the left-hand side of Eq. (19) and the 
value of ey as R, € range arbitrarily. It can be shown 
that ej. is a maximum at the center of the circle and 
further that 
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(€12) maz = byr* (22) 


where 7 is, as above, the radius of the energy circle for 
é2 = 0. It follows that e,. can be represented geome- 
trically as a paraboloidal surface horizontal sections 
of which are the loci of the vector R, €«; see Fig. 1, in 
which the subscript ‘‘em’’ refers to values appropriate 
to (€12)mar. Clearly, Fig. 1 can be interpreted also 
as the influence on the critical energy circle (e12. = 0) of 
variation of the damping ratio or, more precisely, of the 
direct aerodynamic damping coefficient* do». 


Ultimate Stability 


For ej). = 0, the energy circle defines physically pos- 
sible configurations which permit flutter—it represents 
a necessary condition for flutter. Whether flutter 
actually occurs and if so, to which point of the circle it 
corresponds, is determined by system parameters in- 
volving generalized masses and structural stiffness or 
modal frequencies, which are absent from Eq. (19). 
In order to determine the necessary and sufficient con- 
figuration, one must use a variational principle of least 
Action which normally means solving the Lagrangian 
equations of motion for the eigenvalues and then seek- 
ing the eigenmodes. (The setting up of such auxiliary 
conditions using the equations of motion, for example, 
can lead to a complete solution for the eigenvalues as 
demonstrated by Pulos*‘ for a binary system.) Lack- 
ing this knowledge, we note that Eqs. (19) and (21), 
and Fig. 1 indicate that there are some configurations 
represented as lying outside the critical energy circle 
for which the Action is such that the system can never 
absorb energy and therefore never flutter, crrespective 
of the mass and structural stiffness characteristics. It 
follows that if the critical energy circle can be made to 
vanish by manipulation of the significant aerodynamic 
couplings, then stability is always assured. We can 
then define a criterion of ‘ultimate stability’ which 
states that if 


ro” <0 (23) 
flutter instability is impossible. 


In view of Eq. (22) and the fact that the absolute 
values of 


= = 2buro (24) 


where X = Roose, Y = KR sin «, it is plausible to at- 
tempt an interpretation of the significance of the size of 
the critical energy circle. It should first be realized 
that the paraboloidal surface of stability (or instability) 
strictly has no physical validity except at ej. = 0 inas- 
much as the expression for the energy exchange ej. is 
derived from considerations of sinusoidal motion. Let 
us suppose, however, that this surface is an adequate 


* One would expect Eqs. (19) and (20) to show cyclic behavior 
with respect to the subscripts 1 and 2, so that b22 is not a preferred 
quantity. The cyclic nature of the results is achieved by recast- 
ing Eq. (19) by placing Re = R. Succinctly, as formulated the 
analysis deals with R = Q;/Q:2 and not its inverse Q2/Q:, which 
corresponds to a definition of = = dy,/bee. 


ENERGY BALANCE 707 


approximation of the true one which results from a 
correct representation of the motion as exponentially 
damped (positively or negatively). Then any non- 
zero value é,. can be regarded as the quantity of energy 
which must artificially be added to (or subtracted 
from) the system to ensure satisfaction of the equation 
of energy balance and so a critical flutter motion. 
Alternatively, the maximum value of é;:. can be thought 
of as the maximum available energy which the system 
is capable of absorbing. This suggests that the nature 
of any flutter occurring for ¢:, = 0 should be related to 
the quantity (12) mar in such a way that the flutter can 
be described as ‘‘violent”’ if the configuration is poten- 
tially capable of large energy gains and ‘“‘mild”’ if small. 
Through Eq. (24), this means that a large energy circle 
radius 7) implies violent flutter whereas a small r) im- 
plies mild flutter. Thus, the question of the degree 
of violence of flutter is transferred from the usual con- 
siderations of decay rate of amplitude or the rate of 
change of this rate with airspeed say, to the consider- 
ation of the capability of a neutrally stable configu- 
ration of absorbing energy. On the basis of our sup- 
positions, this can be stated in terms of the radius of the 
energy circle so that any parametric change which 
increases/decreases 7, can be described as destabilizing / 
stabilizing. 

We look more closely now at the assumption con- 
cerning the validity of the paraboloidal surface for 


nonzero é). If Eq. (15) is replaced by the motion 


= [Qje" sin (kt + = 1,2,.. 


the equation of energy balance (14), which will now in- 

volve JT and V, yields for the binary system n = 2, 

with R, = 1, Ri = R, = 0, q = A =] + 2d7, Ar = 

AE = [T+ = —[(e* — 1)bukQ:?/2d] X 
[A R? cos? + R?® sin? « — AR? sin 2e + 
2AxR cos « — 2(9E + Ax)Rsin e + Abe] (26) 


.,n (25) 


which is the polar equation of a central conic traced by 
the radius vector R, «in the X = Rceose, Y = Rsin e 


plane. The ellipse center is given by [compare Eq. 
(20)] 

X. = —x + &/(1 + 

Y, = &(1 + 2A*)/(1 + A?) (27) 
the inclination of the major axis to OX is a = —(1/2) 


tan~' (1/A), while the square of the major and minor 
axes are given by 


(1+ + — 1) + 
(1 + 2A?) [ro? — A*E2/(1 + A2)]} (28) 


where T = AE/2rkQ,*; see Fig. 2. The total energy 
exchange can be maximized by considering OF /OX = 
or/oY = or/oaA = 0, the first two conditions of which 
lead to the location in the X, Y plane of Inez; thus 
Xmaz, Ymar are found from Eq. (27) by letting A = 
Amar. The third condition written as 


= 0 (29) 
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Energy balance ellipses for super- and subcritical motion; 
\ = damping ratio. 


Fic. 2. 


reduces ultimately to 


[(1 + + 204)/(1 + 22)? + 
(1 — — + 
Vn/& — (1 + = (30) 


It is realistic now to limit \ so that \ is small compared 
with unity, so that 


(31) 
+ 8(md)? 


and those quantities, x, &, 7) dependent on the flutter 
coefficients may assume their sinusoidal values for which 
rigorously \ = 0. By Eq. (30), it follows that 


Amar = (1/6m)[—1 + V1 + /8)?] ~ 
= + V1 + 

Xmar = + 

(32) 
Expression (28) for the ellipse axes reduces simply to 


(1 + A)[—T/bu(1 + + (33) 


Eqs. (32) form the basis of the claim that qualitatively 
the representation of Fig. 1 (based on \ = 0) and the 
interpretation placed on the significance of the magni- 
tude of the energy circle radius 7 is an adequate one. 
But we may state further that, for small \, the magni- 
tude of the maximum energy capable of absorption 
is a simple function of r) defined by the second of Eqs. 
(32); and, in fact, there corresponds to this condition 
a maximum (positive) damping ratio \,, which the sys- 
tem can attain which is related to 7) in the manner of the 
first of Eqs. (82). We choose then to define the nature 
of a system’s critical flutter oscillation in terms of its 
energy absorption potentiality. 

We conclude then that if in practice it is not possible 
to satisfy the criterion (23), it is at least desirable to 
minimize the energy circle radius by adjustment of the 
parameters involved in x, &, ¢. 


} maz 
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Higher Order Systems and Stability 


For systems with many degrees of freedom, n > 2, 
Eq. (18) does not prove as tractable as for the binary 
case. It is often no restriction to limit considerations 
to a single pair of freedoms as an approximate repre- 
sentation of the fluttering system since such degenerate 
binary forms contain all the essential characteristics 
of the more general system. It is important, of course, 
to select the appropriate parent binary form and this 
is not always an easy choice for today’s configurations, 
Then again it is increasingly the case that a binary 
system is an inadequate approximation for modern 
configurations. It is of interest, therefore, to look 
briefly at the general Eq. (18), and it proves fortuitous 
that this can be rewritten quite rigorously in a form 
involving the m(m — 1)/2 principal binary constituents. 
We define e, as the expression—see, for example, Eq. 
(19)—for energy for the binary constituent formed 
from the system whose two displacements are g; and 
gx, all other g being zero. Then Eq. (18) can be written, 
without approximation, 


n 


e= C ik + (n = 2) b),R;? = 0 (34) 
j=1 k>j j=l 


j= 


in which for all damping coefficients b;; > 0, which is 
the usual case, the single summation is always positive. 
This equation simply states that the total energy bal- 
ance is the sum of the energies of the constituent binary 
systems and the individual modal dissipation energies. 
For simplicity, consider the ternary system, ” = 3: 


= C12 + C13 + C23 + by + + (35) 


Suppose flight conditions are such that all three binary 
systems are unstable—i.e., each e;, > 0; then clearly 
€123 > O and the ternary system is also unstable. The 
same result applies to the general Eq. (34). On the 
other hand, suppose all binary systems are stable— 
i.e., all e;, < 0; the higher order systems are not neces- 
sarily stable and provided some or all of the direct 
(positive) damping coefficients are sufficiently /arge they 
may be unstable. Other corollaries are simply deduced, 
Viz. : 

(a) If all binary systems are critical (e;, = 0), all 
higher order systems are stable/unstable provided all 
dampings b;, S 0. 

(b) If one binary system is critical (say, é12 = 0) and 
both the other binaries are unstable (e)3, ¢2; > 0) for the 
same flight conditions, then for these conditions of 
flight the ternary system is necessarily unstable. In 
other words, the addition of the third degree of freedom 
is destabilizing. However, the converse is not neces- 
sarily true. 

Eq. (34) thus offers a guide to the relevance of higher 
modes. Béatrix!* has suggested a method directed to 


the same end based on a study of the way in which 
constituent modal energies change with airspeed. See 
also Van de Vooren.” 
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Structural Damping 


Structural damping is variously represented.” If a 
viscous type damping force is an adequate representa- 
tion, then these forces appear in Eqs. (3) or (4) simply 
as added components of the damping matrix. Thus, B 
is repleced by B + D where D is a diagonal matrix of 
modal damping coefficients; b;; + d;,; dj, = 0, 
j # k. This substitution carries right through the 
analysis of the energy balance. The more usual “‘hys- 
teretic’’ representation of structural damping—see 
Eq. (2)—-stems from modal damping forces of the form™ 
or dg \q Energy dissipation per cycle of 
oscillation arising from such forces is then, respectively, 


a= 


(36) 
or = 4,0; 
j= 
which differ only in a constant multiplier. Eqs. (16), 


(18), and (19) then possess a right-hand side of the 
form d;,Q;? cos (€; — cos Ajj, and (dyR? + 
d»)/k, respectively, where the constant factor of Eq. 
(36) has been absorbed in d,;. It is possible then to 
absorb the effects of structural damping simply by 
defining new damping coefficients 


by, = bj; d;;/k (37) 


which influence system stability through the parameters 
of Eq. (21) for the energy circle radius. By consider- 
ing O7?/Obe. and Or?/0b, the variation of 7 due to 
changes in by), be. is given by 


aro? = + + (38) 


Assuming, as usual, § > 0, bn > 0, it follows that for 
m2 = 0 increases (decreases) of either or both by and de» 
cause reductions (augmentations) in 7? and, there- 
fore, are stabilizing (destabilizing). [It is of interest 
to compare this result with that deduced from Eq. (35) 
for ternary flutter—there is no paradox since the sta- 
bility determining quantities are essentially x, & of 
Eq. (20), not the direct damping coefficients.] For 
n°? < 0—i.e., for an ultimately stable system changes 
in dampings can cause stabilization or destabilization 
according as 279°/f is S 1. Either possibility exists 
also when 6),, 5b. have opposite signs—see also the 
results of Salaiin." 


Degenerate Binary Systems 


In the unusual event that a system has § = 0—i.e., 
the aerodynamic stiffness matrix is symmetric, the 
energy circle radius becomes 7,2 = x* — £, so involving 
only damping coefficients. Stability questions are 
then related only to the dissipation function F, as pre- 
viously noted. This can be shown in a direct way by 
noting that for m = 2, F is a positive definite form if all 
discriminants Bj | > 0; i2., 


by > O and — + by)/2]? > 0 


which is x?) > 0 
— by (39) 
or, finally, ro” <0 


It is more often the case that x? — ¢ is small and 
negative and £* is so large as to dominate the magnitude 
of ro” ‘In this case, a reasonable approximation for the 
purposes of assessing the system is to place x? — ¢ = 0 
so that the skew-symmetric stiffness components are a 
direct indication of system stability. Mere inspection 
of these flutter coefficients will suggest which binary 
couplings are significant. The neglect** of the damp- 
ing terms leads to a conservative approximation. 


Aerodynamic Balancing 


The application of the principles above to the pre- 
vention or elimination of flutter can be thought of as 
aerodynamic balancing. Early techniques of flutter 
prevention often consisted essentially of dynamic bal- 
ancing—i.e., the elimination or optimization of inertial 
couplings. Here we are primarily concerned with bring- 
ing within acceptable bounds the aerodynamic cou- 
plings. Mass and structural stiffness distribution are 
not explicitly involved, but their implicit role appears 
in the dependence of the aerodynamic coefficients bj, 
c,, on mode shapes as well as their influence on the re- 
duced frequency & at flutter. 

In subsequent sections, we proceed to examine some 
examples of supersonic flutter in order both to gain an 
insight into the essential mechanism and to discover 
in what way useful application of the ultimate stability 
concept may result. 


General Supersonic Flutter 


Consider the regime’ for which piston theory oscilla- 
tory airforces are appropriate. The generalized forces 
can be put in the nondimensional form 


Qik = Ce + = 


in which f; = /;(y, 7) define nondimensional modes as 
functions of the nondimensional coordinates y (stream- 
wise; ¥y = x/2d) and n (spanwise; 7 = y/s) and the 
prime means 0/Oy. IT is a streamwise profile param- 
eter, such that 


= (1/2)(vy + (41) 


when the binomial expression'* for the pressure coeffi- 
cient is truncated to exclude downwash terms of degree 
> 3. In Eq. (41), 20(y) defines the profile shape. It 
follows* immediately from Eq. (40) that the dissipation 
matrix [b,] is symmetric and that for f,’ = 0; ice., 
for any pure bending mode, c,, = 0 so that there exist 

* The only restriction on profile shape implicit in the derivation 


of the form (40) is the requirement of profile symmetry about the 
mean chord-line. ‘‘Fore and aft’? symmetry is not necessary. 
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Fic. 3. Energy circles for binary wing bending-torsion flutter. 
(a) Supersonic flutter compared with incompressible for midchord 
pitch axis of two-dimensional section. (b) Effect of aspect ratio 
on flutter of rigid flat plate at MW = 2. 


likely large skew-symmetric components ¢;, — In- 
terpreting this in binary terms, we find that 


X = = — = Dex (42) 


For a surface of (near) uniform mass distribution, the 
orthogonal modal coupling integral f, is (small) 
zero so that for no thickness effects, f = 0, bj, = x = 0. 
It can be expected then that the skew-symmetric 
effects embodied in £ will dominate the magnitude of 
the energy circle radius for any binary coupling in- 
volving a bending type mode, in such a case. This 
result is closely true for more general modal configu- 
rations. 


Binary Bending-Torsion Wing Flutter 


If the above expressions are developed for the simple 
two-dimensional, infinite aspect ratio wing with a bend- 
ing mode f; = 1, and torsion mode fz = x/2b — a, where 
x = 2ba defines the pitch axis position aft of the leading 
edge we find 


x = 1/2 —a — (y + 1)MA,/16b? 
= 1/4k 
1/3 —a +a? — (y + 1)MA,(1/2 — a)/86? 


(43) 
sothat ¢ — x? = 1/12 — [(y + 1)MA,/16b?]* (44) 


and 
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By inequality (39), the dissipation F is positive 
definite, from Eq. (44), provided (y + 1)MA,/b? < 
8/+/3, which for a doubly symmetric wedge profile in 
air reduces to Mr < 5/3+/3 ~ 1. It is clear then that 
for practical profiles F is a positive definite form, ¢ — x? 
is small compared with £? for the small values of k ex- 
perienced in this type of flutter, and the energy circle 
radius is predominantly determined by £. For the 
flat plate profile, specifically, and rigorously, 


r% = (1/2)V1/4k? — 1/3 (46) 


a result independent of both Mach Number and pitching 
axis, which vanishes for k = k,, > +/3/2 when flutter 
is impossible, the system having attained ultimate sta- 
bility. For the general profile, it transpires that the 
influence of thickness on the energy circle is very small 
for a prescribed k. It is always mildly destabilizing 
when & is invariant with change in the thickness param- 
eter WA,,/16b?. However, in practice, changing thick- 
ness normally results in a change in & the nature of 
which is not revealed by the equation of energy balance 
approach. The qualitative effect can be estimated by 
looking at the change in the modal aerodynamic stiff- 
nesses which will influence the flutter frequency and 
therefore k. Thus, c, = 0 and cx = [4(1 — 2a) — 
(y + 1)Mr]/4M for the double wedge section which 
decreases with increasing \/7. Accordingly, k can be 
expected also to decrease so leading to increased 1% 
and destabilization. For a pitch axis near midchord, 
a— 1/2, the effect may well be significant."* 

It is of some interest to compare the results (43) and 
(46) with those arising from another quasi-steady 
theory,'® viz., two-dimensional incompressible flow. 
By arbitrarily placing the Theodorsen function C(k) = 
1/2, 1 we derive two approximate formulations of the 
critical energy circle. We find in both cases that like 
the piston theory result (46) the radius is independent 
of pitch axis position whereas the coordinates of the 
center are given, respectively, by [(1/2 — a)/2, 1/4k], 
[1/4 — a)/2, 1/4k]—compare Eq. (43) for zero 
thickness. The energy circles are graphed in Fig. 3(a) 
for the pitch axis at midchord and k = 0.3, where, in 
addition, is given that for the exact linearized theory 
for M = 0. We note that for C(k) = 1/2 the energy 
circle passes through the origin since § = 0, where R = 0 
and ¢ is indeterminate. (This result is independent of 
k for a = 1/2.) This is interpretable as a zero air- 
speed natural oscillation of a single mode. 

Fig. 3(b) illustrates the effect of aspect ratio on the 
energy circle for a rectangular flat plate fluttering in 
rigid translation and pitching about the midchord. It 
is assumed that for AX = © piston theory is acceptable 
at MM = 2. The finite M circles are based on correc- 
tions to the R = © case defined by a linearized theory” 
which is considered valid for R > 1/+>/ M? — 1 = 0.58. 
Decreasing A is seen to be destabilizing in this case. 
In Fig. 4 are given further results for pitching about 
the leading edge (a = 0); Fig. 4(a) depicts the effect 
of Mach Number for R = ©, while Fig. 4(b) shows the 
influence of R for M = 1.4, k = 0.3. In both casés, 
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exact linearized theory isused. It appears that M vari- 
ation has little influence on the energy circle radius 
whereas decreasing A is strongly stabilizing. 


Binary Wing-Aileron Flutter 


Let 2bh define the chordwise location of the aileron 
hingeline aft of the wing leading edge. Consider first 
wing translation-aileron rotation coupling according to 
piston theory in which now an additional downwash 
term is allowed in the expression for the pressure co- 
efficient so that effects of mean aileron displacement are 
included. It is found that, for a double-wedge profile, 


(1 — h)2?/C 

x/2k(1 — h) 

2x(1 — h)/3 

= —(y+ 1)Mr+ B/2 

(y + 1)M?(r? + By?) 

(4+ B)/(2 + A) (47) 


by which the energy circle is defined (8) = initial aileron 
deflection). Stability boundaries defined by ~™? = 0 
are shown in Fig. 5 as a function of aileron chord (1 — h) 
and reduced frequency k, for the zero thickness case, 
A= B= 0. Fork 3/4, for example, it is impos- 
sible to achieve ultimate stability. 

When thickness and initial aileron deflection are 
allowed with k invariance, it is possible to deduce that 
the former by itself is stabilizing, the latter mildly de- 
stabilizing. Together the overall effect is destabilizing 
only when { is large compared with r. 

Consider now the case of wing pitch-aileron rotation 
coupling for the double wedge profile. Analogous to 
Eq. (47), we then have 


= (1 — h)*(h —a + 2(1 — h)/3)/D 

= (1 — h)(h — a)/2kD 

= 2(1 — h)*/3D (48) 
C(a? — a/2 + 1/12) + (1 — 2a)/2 

= (4 + B)/(2 + A) 

and the vanishing of the energy circle, 7)? = 0, leads to 


a functional relation between the hinge-line location 
hand the pitch axis a, thus: 


| 


ll 


3(1 — Aja? — 211 — h — h*)a — 
(1 +hA+ h?)/3 — -—a+1/6]D + 
3(a — h)?/4k2(1 — hh) = 0 (49) 


Graphical solutions of this condition for the flat plate 
A = B = O are given in Figs. 6 and 7. The former 
applies to a quarter-chord aileron, h = 3/4 and the 
shaded area indicates ultimate stability, 7? < 0. The 
contours r? ~ (0 reveal the existence of a worst siting of 
the pitch axis, viz.,@ = 1/4 where for a given k, the 
energy circle radius is a maximum. If increasing k is 
regarded as being due to increasing a modal stiffness 
(the aileron circuit, say), then the contours suggest that 
in the region a = 1/4 such a measure is not effectively 
stabilizing. Rearward positions of the pitch axis appear 
more favorable. Fig. 7 summarizes the critical condi- 
tion for0 < h 1,0 <a 2for a variety of k, reflect- 


PISTON 
THEORY 


(a) 


Fic. 4. Energy circles for binary wing bending-torsion flutter. 
(a) Effect of Mach Number on two-dimensional section pitching 
at leading edge. (b) Effect of aspect ratio at M=1.4 for rigid 
flat plate pitching at leading edge. 


ing a general behavior similar to that for h = 3/4, 
Fig. 6. For h = a, the system possesses ultimate 
stability for all k. 


Flutter of All-Moving Surfaces 


The configuration used here consists of a finite wing 
symmetrically tapered about the locus of midchord sta- 
tions, having a taper ratio ¢ = tip chord/root chord. 
Two modes are considered, the first being rigid body 
pitch about an axis 2ba aft of the root leading edge, the 
second being a bending mode represented by the non- 
dimensional modal function f = — 4?)/(1 — 
The value 4 of the spanwise coordinate 7 thus defines 
the spanwise position of the bending nodal line. For 
the delta wing, ¢ = 0, the use of piston strip theory will 
yield correct values*! of many of the generalized forces 
for the modes considered, while those remaining will be 
closely correct. For ¢ > 0, the approximation is less 
good because of the parallel tip effect, but for the pur- 
pose of illustrating the principle of ultimate stability 
this error is neglected. 

Fig. 8 depicts the critical stability surface, 7? = 0 
in terms of pitch axis location a, bending nodal line 
position 7, and reduced frequency k based on the root 
chord for the delta wing ¢ = 0. The surface is sym- 
metrical about the 4:k plane at a = 1/2 which for in- 
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Fic. 5. 


(1 — h) and reduced frequency. 


PITCH AXIS POSITION 


Fic. 6. Binary wing torsion-aileron supersonic flutter showing 
contours of the energy circle radius squared as function of pitch 
axis position and reduced frequency for hinge position of 3/4 
chord. 


DELTA WING 
t=O 


MS 


Fic. 8. Supersonic binary flutter of all-moving delta showing 
stability surface as function of nodal line and pitch axis positions. 
Ultimate stability exists above the surface. 


Binary wing bending-aileron supersonic flutter show- 
ing ultimate stability boundaries as function of aileron* chord 
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HINGE-LINE POSITION (h) 


PITCH AXIS POSITION (a) 


Fic. 7. Binary wing torsion-aileron supersonic flutter showing 
limiting reduced frequencies for ultimate stability as function of 
pitch axis and hinge-line positions. 


CRITICAL REDUCED 
FREQUENCY k,, (7 =O) 


k 
| -TTCRITICAL NODAL 
LINE POSITION fe 


| | 

J 
| = 
© TAPER RATIO (t=b,/b,) 


Fic. 9. Flutter of all-moving surfaces showing variation of 
critical nodal line position and reduced frequency for 7,= 0 with 
taper ratio. 
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(,) 


Fic. 10. Flutter of carried all-moving rectangular control 
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REDUCED FREQUENCY (k) 


showing limiting reduced frequency as function of pitch axis 
position a and contro] span. 
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board locations ot the bending nodal line is clearly the 
worst situation for the pitch axis. Ultimate stability 
exists above the surface. It is of interest to note that 
there exists an » = 7., for which ultimate stability is 
assured irrespective of k. Its value in this case is 
te = 1/6 and is the optimum nodal line location— 
the system will never flutter. The effect of taper ratio 
on the value of 7,, is shown in Fig. 9 where also is plotted 
for 7 = 0 the value of k,, which always occurs at the 
worst pitch axis position (a = 1/2). Increasing taper 
ratio is seen, in this case, to be mildly stabilizing. 

When the all-moving control is carried by a wing 
surface, the binary flutter is conveniently formulated 
in terms of the same bending mode and rigid control 
pitch outboard of the spanwise station m) which defines 
the control root chord position. For the case of 7 = 
)—i.e., the nodal line of the bending mode is at the 
wing root, and a zero thickness profile, the critical 
stability boundaries are shown in Figs. 10 and 11. The 
former relates to a rectangular planform, the latter to 
adelta. Ultimate stability exists above the boundaries. 
In both cases, it is clear that from the ultimate flutter 
stability viewpoint the carried control is superior to the 
all-moving wing only for a limited range of extreme 
outboard positions of the carried control. For the 
delta wing, for example, ultimate stability can never 
be achieved for a control span one quarter of the wing- 
span—i.e., for 7 = 3/4. 


Panel Flutter 


Let us consider the two-dimensional panel subject to 
piston theory type forces for which we shall assume 
zero thickness. Further, we use a modal representa- 
tion of the displacements so that for simply-supported 
edge conditions (normal to airstream) the mth mode 
is simply 

fn = sinmr0, m2 1 (50) 
and we make no distinction between a ‘“‘membrane”’ 


and a “‘plate’’ panel. The streamwise coordinate is 
here x = 266, b being the panel semilength, so that 


VA 


Ww 


REDUCED FREQUENCY (k) 


| 
PITCH AXIS LOCATION (a) 


_Fic. 11. Flutter of carried all-moving delta control showing 
limiting reduced frequency as function of pitch axis position a 
and control span. 


CRITICAL REDUCED FREQUENCY (k¢,) 


\ 


MODES m,n =1,2,3 


Fic. 12. Two-dimensional panel flutter showing limiting 
reduced frequency for binary couplings of modes f,, = sin m6 
and f, = sin 


the panel is defined between 6 = 0, 1. Using Eq. 
(40), it is readily shown that the aerodynamic damping 
coefficients are (omitting the common multiplier 2/1/) 


Dum = 1/2; Onn = 0, all m, n (51) 


The damping matrix is thus diagonal. The aerody- 
namic stiffness coefficients are found to be 


Cmn = (0; m # n, both odd or even 
2mn/(m? — n?);m+nodd (52) 
Cnn = 10 


The skew-symmetric stiffness components, therefore, 


are 


Cmn® = (Can — Cam)/2 = 2mn/(m? — n*) (53) 


For any binary coupling, the energy circle radius be- 
comes 


ro? = 4[mn/(m? — n*)|?/k? — 1 (54) 


and the critical value of reduced frequency correspond- 
ing to ultimate stability is derived from 7? = 0 as 


= +2mn/(m? — n?) and m? = (k/k)? — (55) 


Eqs. (51) and (52) indicate that the results (52)—(55) 
apply only for m even, m odd, or vice versa. In all 
other cases, the modes are completely uncoupled, and 
flutter considerations do not arise. 

Under prescribed airstream conditions or for pre- 
scribed k, Eq. (55) shows k,, to be a direct indication of 
the energy circle radius and, therefore, of the signifi- 
cance of coupling of any two modes. In Fig. 12 are 
plotted* values of k,, for some of the lower modes, 


*In Fig. 12, the continuous curves are for convenience of 
visualization only; m and n are, of course, integers. 
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and it is at once evident that whereas for increasing n 
the energy circle radius rapidly decreases (toward zero), 
for increasing m the radius increases rapidly for low n. 
This suggests that such higher modes are of vital signifi- 
cance in a modal-type analysis but it is not clear how 
this significance is manifest; and it is of interest that 
the analysis so far takes no cognizance of any difference 
between the behavior of a membrane and a plate. 
The discrimination of the two structures arises solely 
in the disparity of their respective actual flutter reduced 
frequencies k in Eq. (55) due to differing modal fre- 
quencies for the same mode. This is most clearly seen 
by an assessment of the flutter reduced frequency k, by 
means of a conventional solution of the binary flutter 
determinant. The separated characteristic equation 
yields finally, for any two modes m, n: 


= (Rm? + (56) 
where k,, and k, are reduced modal frequencies, and 


where yp is the density ratio, the exact definition of which 
is immaterial. Eq. (56) allows discrimination of a 
membrane and plate by use of the k,, and k,, appropriate 
toeach. Eq. (57) shows that for flutter to exist 


ke < ke? (58) 


which is exactly the condition derivable from the cri- 
terion of ultimate stability through Eq. (55). In 
other words, ultimate stability is here synonymous with 
stability in the usual sense. From Eqs. (56) and (57), 
provided Eq. (58) is satisfied, the flutter speed is given 
by 

U, = — Wn”) /w (59) 


Consider now the values of k,, Eq. (56). The modal 
frequencies and, therefore, modal reduced frequencies 
k,, are proportional to m and m?’, respectively, for the 
membrane and plate. The constants of proportion- 
ality depend on membrane and plate mass and stiffness 


[ox] =] 1-1 0.59 —0.28 
0 0.25 —0.30 0.14 

9 0.047 —0.095 

0 O 0.012 0.011 

0 0 —0.015 0.0022 

[Cex] =[O -1.5 1. 


0 0 
0 0 


0 0) — 


together with a symmetric matrix of virtual inertias 
and a diagonal generalized inertia matrix. Our object 
is to assess the relevance of each mode and the likeli- 
hood of flutter. Form first the symmetric damping 
matrix by = (bj + b:;)/2 and, secondly, the skew- 
symmetric stiffness matrix cy? = (¢, — Cx;)/2, which 


0 


—().13 Ist overtone elastic mode (4) 
0.21 | 2nd overtone elastic mode (5) 
8 -—1.6 9.9 
3 0.36 —4.2 
2.0 
.05 —0.088 —0.57 
05 0.066 0.57 


properties. Let these be so chosen that for both plate 
and membrane k, assumes the value of &,,, for binary 
coupling of the lowest modes—i.e., k, = k,- = 4/3 for 
both structures for m = 1,n = 2. It is simple to show 
then that for the plate k,(m, n) > k.,(m, n) for all other 
m, n; while for the membrane k,(m, n) < k.,(m, m) for 
all couplings m, n = m + | but exceeds k,,(m, m) for all 
remaining binary couplings. Interpreting this in terms 
of Fig. 12, the ‘‘Surface’’ formed by k, for the plate 
lies everywhere above the ‘‘Surface’’ of k,, shown in 
that figure, excepting at m = 1, m = 2. In other 
words, all the higher mode binaries are stable. For 
the membrane this is not so, and for all binary cou- 
plings corresponding to the line AB of Fig. 12, k, < k,,, 
and all such binaries are unstable at some finite air- 
speed; elsewhere all binaries are stable for the mem- 
brane. The two structures are thus basically different, 
and it appears intuitively that whereas for the plate 
the higher modes become increasingly less significant, 
for the membrane a// the higher modes attain impor- 
tance in a Rayleigh-Ritz type analysis.'* Unfortu- 
nately, the results above do not allow a rigorous demon- 
stration of this. 


Chordwise Flexible Surfaces 


The flutter of an extremely small aspect ratio delta 
wing provides an examplet of the limitations of the 
energy balance approach when dealing with multi- 
degree of freedom systems as well as highlighting the 
dangers of disregarding structural parameters. The 
wing is regarded as possessing no spanwise deformation, 
the vibration modes being polynomial functions only 
of the streamwise coordinate. Calculated by an iter- 
ation process, the first three elastic modes are con- 
sidered and shown sketched in Fig. 13. In addition, the 
two rigid body modes of translation and pitch are in- 
cluded in the analysis. The modes, designated 1, 2, 

., 5 have modal frequency parameters which can be 
represented, respectively, as 0, 0, 1, 3, 8 approximately. 
When the flutter coefficients are calculated on the basis 
of slender body theory,” we find 


0.98 | rigid translation (1) 
—().49 | rigid pitch (2) 
0.31 fundamental elastic mode (3) 


is easily done mentally. We notice that the former are, 
in general, quite small compared with the corresponding 
latter elements. which can, therefore, be expected to 


t Taken from work currently sponsored by the AFOSR under 
contract No. AF49(638)-219. 
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dominate the magnitude of the energy circle radius for 
the various binary couplings between the five modes. 
It is clear too that the skew-symmetric components for 
the rigid-body modes (the first and second rows and 
columns) are significantly large. On this basis then, 
we are led to conclude that both rigid body modes can 
be expected to play a leading role in flutter of the 
coupled system. However, this in fact is not so, 
simply because all subdiagonal elements of all coeffi- 
cients a (inertia), b, c, as well as the virtual mass co- 
efficients are zero in columns 1 and 2—.e., for the two 
rigid body modes. This means that in the flutter de- 
terminant all the corresponding 6;,—see Eq. (1)— 
are zero so that these two modes are completely un- 
coupled from the elastic modes; they play no part in 
the flutter of the coupled system. The energy balance 
concept thus makes no explicit allowance for such 
partly-triangular flutter determinants; it essentially 
provides a necessary condition for flutter. 

Deleting then the first two rows and columns of the 
coefficient matrices, we assess the remaining binary 
couplings, 3-4, 4-5, 3-5. A little mental arithmetic 
shows that the energy circle radii of the first two are 
of the same order whereas that of the third is somewhat 
larger, which suggests that mode 5 has an important 
role and should be included in a conventional flutter 
analysis. Such an analysis, however, reveals that the 
3-5 and 4-5 binaries are always stable, that the 5-4 
binary produces critical flutter conditions which, as 
may now be expected, differ slightly and only quanti- 
tatively from the critical flutter conditions for the 3-4-5 
ternary coupled system. We conclude that although 
mode 5 is potentially capable of a flutter coupling with 
mode 3, the structural parameters are inappropriate to 
allow it—the ratio of modal stiffnesses, for example, is 
(8/1)?. 

This is not strictly a fair test inasmuch as the energy 
circle concept rightly is restricted to binary systems. 
More rigorously, Eq. (35) should be consulted in con- 
junction with the known results of the binary flutter 
analyses. Thus, for flight conditions which make the 
3-4 binary critical, e;; = 0 and 


eus = — + 0.05Rs? + 0.01R2 + 0.2R3? 
(60) 


since 4-5, 3-5 binaries are always stable. It is impos- 
sible to decide whether e3;; 2 0; i.e., whether the in- 
clusion of mode 5 stabilizes or destabilizes the critical 
3-4 binary. 


Some Transient Effects 


We consider a binary system and make the assump- 
tion, ab initio, that the introduction of time-dependence 
of the flutter coefficients need not cause deviations of 
the type of system motion—i.e., the motion (15) is 
valid. Substituting Eq. (15) into Eq. (13) and spe- 
Cializing to n = 2, we find 


ENERGY BALANCE 715 


MODE | 


Fic. 13. Small aspect ratio chordwise-flexible delta wing modes, 


2n/« 2a /K 
0 
R? + 2(x cose —tsine)R+¢} (61) 


so that structural parameters now enter the equation 
of energy balance. To give some reality to the dis- 
cussion, Eq. (61) is applied to a system subject to aero- 
dynamic heating effects which are assumed to be mani- 
fest only by a time rate of change of the structural stiff- 
ness coefficients. Thus, we place A = 0 and the ele- 


ments of ¢ equal to 


Cix = 0,7 x k 
= psye? (62) 


where o;;, sj; are aerodynamic and structural stiffness, 

respectively. Eq. (61) reduces to 

—(T + V)?"* = aR*® + BR* cos 2e + rR? sin 2e + 
6bReose + @Rsine +1 = 0 (63) 


where 


= 2rkbo» P02 

+ (p/2k)?] 

= —(p/2k)B/2 

2rk (dye + be) 

— Psufl + 2(p/2k)/[4 + (p/2k)?]} 
P = le 2x(p/2k) __ 1]/4 


WR 


II 


which is the polar equation of a conic in the X¥ = R cos 
e, Y = Rsin e plane. Its nature is more easily seen 
if the restriction is made that 27(p/2k) is small—i.e., 
if time variations of stiffness are not great, an assump- 


tion which accords with usual circumstance. Then 
= 1 (p/2k) 5022/2] 
== 
(64) 


= (die + 
= 29 (C12 — C2) 


m(2kdu (p/2k)su/2] 


II 


and Eq. (63) defines kinematically admissible flutter 
configurations by means of a second energy circle the 
characteristics of which differ from that of the conven- 
tional steady-state case only in relation to effective 
damping coefficients. 


Thus, by defining 


‘ 
| 
| on 
| | | 
| 
i 


= bj; — 7 = 1, 2 (65) 


the same energy circle is valid. Time-dependence 
thus appears as an effective damping which is positive 
or negative according as p S (—i.e., as the stiffness 
decreases or increases. An interesting result. 

A simple application is to the two-dimensional binary 
flutter of a zero thickness wing for which a rigid pitch 
mode couples with a parabolic chordwise bending mode 
according to piston theory. The limiting condition 
r? = 0 leads to 


(1 — 2a)? + (1/k)? = 36[2(1/3 — a + a?) — 
Sup/8k*] [1/3 — sexp/8k*] (66) 


for which it 7s possible to find real values of p. Fora = 
1/2 for example, p is real provided 


(su $o9/2)? + (1/k? 2) Pg 0 (67) 


Finally, we consider as a second example the two- 
dimensional wing bending-torsion flutter using zero 
thickness piston theory again. For a mid-chord pitch 
axis, the effective dampings are 


(budets = (2/M) — (p/8k*)(uM/3)ko? 
(bee)ett = (2/3M) — (p/8k*)(uM/3)ko/20 


where ky and WV define the stiffnesses 
= (1/3) C2 = Ve, /20 


The condition for the vanishing of the energy circle, 
ro? = 0, results in 


(4/M — puMky?/4k2)(4/M — = 
12/(Mk)2 (68) 


Fig. 14 depicts a solution for the value of the modal 
stiffness ratio WV = 5/3, and illustrates the stabilizing 
influence of a negative rate of change with time of 
stiffness. 


Conclusion 


The equation of energy balance is but a partial in- 
tegration of the equations of motion and as such is in- 
capable of supplying the eigenvalues or the necessary 
and sufficient conditions for the existence of flutter. 
It does offer a means of deriving a necessary condition 
and as such provides an insight into the basic mecha- 
nism of flutter and an understanding of the significant 
characteristics of fluttering systems. For binary sys- 
tems, it leads to the idea of ultimate stability and a con- 
cept of stability related to the capability of a neutrally 
stable system of absorbing energy. A criterion results 
which exposes in a specific way the influence of certain 
system parameters, and this requires as data simply 
the flutter coefficients as conventionally derived. 
For higher-order systems there appears to be no such 
simply applicable tool of assessment, although system 
behavior is often then approximated by certain binary 
couplings to an extent adequate for purposes of pre- 
liminary design. 

Finally, such nonclassical effects as time-dependence 
arising from changing mass, altitude, aerodynamic 
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Fic. 14. Effect of time change of structural stiffness on two- 
dimensional supersonic wing bending-torsion flutter showing 
limiting reduced frequency as function of a stiffness variation of 
the form e?*. 
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heating are catered tor; and there seems no reason why 
the influence of other force fields and other work func- 
tions, linear or nonlinear, should not be added to the 
conventional. 
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Some Effects of Sound-Reduction Devices on 
a Turbulent Jet 


GILLES M. CORCOS* 
Douglas Aircraft Company, Inc. 


Summary 


An account is given of an attempt to explain in simple terms 
how ejectors, corrugated nozzles, and multiple nozzles modify 
the turbulent mixing of a jet so as to generate less aerodynamic 
noise. Simple analytic and dimensional considerations lead to 
the conclusion that the noise reduction is due primarily to a de- 
crease in turbulence levels in the mixing region. It is suggested 
that this is accomplished by accelerating the ambient air through 
a longitudinal pressure gradient before mixing it with the jet. 
Thus, a sheath of flowing air decreases the lateral jump in the 
mean axial velocity which is responsible for the turbulence level 
in the mixing region. An experiment at low speed with an ejector 
and a corrugated nozzle is reported. Mean velocity profiles and 
turbulent surveys were obtained. The results are in good agree- 
ment with expectations and provide strong support for the 
hypothesis. 


Introduction 


: hws DESIGNER of a jet exhaust nozzle, attempting 
to minimize jet mixing noise without incurring 
performance penalties—i.e., loss of thrust and an 
increase in drag and weight—has had little help from 
analysis or from past experience. His approach has 
therefore been purely pragmatic. Considerable results 
have been achieved and it is of interest, in the wake of 
this progress, to seek a measure of understanding. 

It is established that the major part of the noise 
generated by jet engines originates in the mixing of the 
high velocity exhaust gases with the air surrounding 
the jet. This mixing is turbulent and characterized by 
large velocity fluctuations and consequently by pres- 
sure fluctuations which in part propagate or radiate 
away from the jet. The magnitude of the pressure 
fluctuations in the jet is related to the magnitude of 
the velocity fluctuations there, and, as is shown in ref- 
erence |, the extent to which these pressure fluctuations 
propagate away depends on their frequency or time 
scale. 

Several years ago, M. J. Lighthill? provided an 
analysis which yields the radiated noise at a point (or 
its stochastic description) in terms of the kinematics of 
a turbulent field such as a jet mixing region causing the 
noise. From an analytic standpoint this would solve 
the problem, were it not for the fact that the turbulent 
field which is assumed given in Lighthill's analysis 
is, in reality, not known. In other terms, our source 
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of noise is identified but not given by that analysis or 
by any other. Since the problem is in fact not yet 
mathematically determinate (see, for instance, ref- 
erence 3), the designer should not wait upon future 
progress in turbulence theory. Is it possible, instead, 
to measure the relevant quantities in the turbulent 
region and relate them to the geometry of the nozzle 
and the flow parameters? Here it seems that the terms 
identified in Lighthill's theory as the source terms can- 
not yet be measured practically’ even if we overlook 
the difficulty of interpreting hot wire signals at the 
relevant mean velocities—i.e., transonically. 

There are alternate identifications of the source (in 
terms of density, for instance) which may lead to 
practical measurements of its strength by optical 
means. Such measurements have not yet been re- 
ported, however, and in the meantime it is natural to 
seek a substitute—i.e., to associate more approximately 
the strength of the source of noise with a quantity which 
is readily measured in the mixing region and on which 
the effects of changes in the geometry of the nozzle can 
be easily traced. Such a substitute is proposed later. 


Successful Devices 


We shall briefly review modifications of exhaust 
nozzles which have been found successful to various 
degrees in that they generate less noise for a given 
thrust and mass flow than plain round nozzles do. 

It has been verified that an ejector placed immediately 
at the exit of a plain nozzle is capable of decreasing the 
radiated noise energy by a factor as high as 5. This 
occurs only when the stream issuing from the ejector 
has reasonably uniform total pressure—t.e., when 
secondary and primary air have been sufficiently mixed. 
“Sufficient’’ mixing requires a long ejector (the ejector 
length should be approximately eight times its diameter). 

Nozzles, the perimeters of which have been deeply 
corrugated so as to form lobes have been found much 
less efficient noise generators than comparable round 
nozzles. Here the gains are again of the order of 7 
db. or a factor of 5. 

Comparable figures have been quoted if the exhaust 
jet is subdivided in a number of parallel round nozzles 
placed a suitable distance apart from each other. 

The gains obtained by ejectors and by corrugated 
nozzles are to a certain extent additive, in that a com- 
bination of a corrugated (or multiple) nozzle and an 
ejector reduce the aerodynamic noise of a jet still 
further. In fact, it has been found that even very 
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Fic. 1. 


short ejectors were effective complements of corrugated 
nozzles.‘ 


The Dependence of Aerodynamic Noise 
on Simple Flow Parameters 


To help us discuss the role of various parameters 
characterizing our turbulent mixing process and the 
resulting noise (such as the lengths associated with the 
nozzle geometry, the mean and the turbulent velocities, 
characteristic frequencies, etc.) we will use the rela- 
tion derived inreference 1, a low Mach Number ap- 
proximation : 


V*p — = (1) 


Eq. (1) gives the noise pressure field p(x, y, z, #4), in 
region II far enough away from the region of turbulent 
mixing, in terms of the local turbulent pressure fluctua- 
tions Pi»; which are found in the region of turbulent 
mixing (region I) and which are closely given by in- 
compressible relations (see Fig. 1). The right-hand 
side of Eq. (1) is the source of the radiated noise, and 
it can be seen that this term is large either when local 
turbulent fluctuations are large or when they are rapid. 
In fact, one should expect from the well-known solu- 
tions of the inhomogeneous wave equation (in the ab- 
sence of boundaries) that the mean square noise pres- 
sure p? at a distance will be proportional to the mean 
square of the local turbulent presssure and to a fourth 
power of a typical turbulent frequency w*. Thus, 


(1/ao*)w**p? ine 


Now pine Will be taken proportional to u’* (where wu’ 
= Vu? is the root mean square of the axial component 
of the fluctuating velocities) as is suggested in ref- 
erence 1. In the present dimensional argument, no 
differentiation is made between the various components 
of the turbulent velocity fluctuations, and it is assumed 
that u’ is sufficiently representative of the two other 
fluctuating components. Thus, we write 


~ 


Therefore, the noise generated can be expected to be 
reduced by reducing turbulent levels. 

The typical frequency, w*, is undoubtedly also re- 
lated to velocities in the turbulent region, but addi- 
tional considerations must be adduced to decide whether 
the turbulent or the mean (say, jet exit) velocities are 
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relevant. The distinction is important since a noise 
suppressing device is not permitted to alter appreciably 
the mean exit velocity (which would change the thrust 
or the mass flow of the engine). If, as Lighthill 
argues, the turbulent velocities determine the fre- 
quencies which are relevant to the radiation of noise; 
i.e., if w* ~ u’/L, where L is a characteristic dimension; 
then the sound level at a distance and, therefore, the 
total sound power output are indeed very sensitive 
functions of the turbulence level. since then p? ~ wu’. 
Thus, it is reasonable to inquire first whether successful 
nozzle modifications result in a general decrease in the 
turbulence levels. 


The Role of Length Parameters 

The most obvious feature of corrugated and multiple 
nozzles is that they modify the length scales of the mix- 
ing process. For the same mass flow and the same 
thrust the ‘“‘wetted perimeter” is increased, the axial 
distance over which a given amount of surrounding 
air is entrained is decreased, eddy sizes are decreased, 
etc. This may be misleading. It can be shown, for 
instance, that gains (reduction in noise output) realized 
by the substitution of multiple nozzles for one nozzle 
have very little to do with this change of scale. 

A simple application of dimensional analysis indicates 
that for geometrically similar jets, the acoustic effi- 
ciency—i.e., the ratio 


Acoustic energy radiated by jet/Unit time 
Kinetic energy flux of the jet 


depends on the diameter of the nozzle at most through 
the effect of changes in the Reynolds Number. We 
should expect that effect to be small over a wide range 
of Reynolds Numbers. Over the range of interest, 
Tyler’s experiments® as well as those of Lassiter and 
Hubbard’ have shown that changes in the Reynolds 
Number do not affect the acoustic efficiency of similar 
nozzles. Thus, if ~ small round jets have the same 
exit velocity and total exit cross-sectional area as one 
large round nozzle, they will generate together exactly 
as much sound power, provided they do not interfere 
with each other; this in spite of the fact that all length 
scales are different. Therefore, it is clearly the inter- 
ference between contiguous nozzles which is respon- 
sible for the reduction in sound power, not the change of 
scale. 

Thus, we have been led to the expectation that a 
decrease in the acoustic efficiency of a jet is closely 
related to a decrease in the turbulence levels within the 
jet. Now assuming this inference to be correct, why 


An ejector accelerates secondary air before mixing it. 
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would, say, an ejector depress turbulence levels? 


It was surmised that what determines the turbulence 
levels is the lateral mean jump in the axial velocity 
between the jet and its surroundings. For instance, 
if a nozzle discharging a jet at a velocity V; were dis- 
placed axially in the opposite direction at the same 
speed, it seems intuitively clear that the turbulence 
levels would practically disappear. Now in an ejector 
the secondary stream is first accelerated as a potential 
flow by a longitudinal pressure gradient, a process which 
causes no turbulence and which decreases the velocity 
of the primary jet relative to its surroundings. The 
subsequent mixing is therefore characterized by a 
lower effective velocity so that turbulence levels might 
be expected to be lower than they are at corresponding 
stations of a free jet. A similar explanation can be 
advanced for multiple or corrugated nozzles. Again, 
the unmixed air surrounded by separate jets will 
acquire an axial velocity which will reduce the relative 
velocity of the jet with respect to its surroundings. 
(See Fig. 2.) 

The foregoing hypotheses were tested in a low-speed 
experiment which will be reported presently. 


Experimental Equipment and Methods 


A one-horsepower vacuum cleaner discharged air 
through a diffuser to a cylindrical settling chamber. 
A smooth Fiberglas contraction section ended the set- 
tling chamber and formed a round nozzle 1.25 in. in 
diameter. In the first part of the experiment an ejector 
14 in. long and 1.9 in. in diameter was secured to the 
end of the nozzle. The ejector was a thin piece of 
aluminum tubing properly flared out at the entrance 
end. 

In the second part of the experiment, a four-lobe 
corrugated nozzle was used. An insert was made out of 
a piece of brass tubing the outside diameter of which 
was precisely equal to the inside diameter of the nozzle. 
One end of the tube was hammered into a four-lobe 
nozzle. The nominal Reynolds Number based on the 
diameter of the basic exit nozzle remained 40,000 for 
all the measurements. Both for the ejector and for 
the corrugated nozzle tests, the boundary layer inside 
the nozzle was made turbulent by tripping it with a 
thin cellophane tape, approximately one inch upstream 
of the nozzle exit or of the start of the corrugations. 

Apart from a few total pressure traverses, the experi- 
ment consisted of hot-wire measurements of the tur- 
bulence levels u’/U (here u’ is the root mean square of 
the fluctuating velocity in the x direction, U the mean 
velocity at the same point) and of U. 

The hot-wire probes were small u’ wires (single- 
wire stretched perpendicularly to the mean flow direc- 
tion). Their length was 0.1 in. and their diameter 
0.0001 in. The measurements were made by the con- 
stant temperature method. The amplifying and aver- 
aging equipment, the bridge, the compensation, and 
the data reduction followed conventional lines. 

Even though the wires were calibrated and the cali- 
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bration checked frequently, all measurements were 
taken in a comparative manner. A standard traverse 
was made immediately after each one or two traverses 
and compared to independent measurements of the 
same quantity. Thus repeatability was assessed and 
changes in the calibration easily detected. On the 
whole, the accuracy was satisfactory. It was found 
that most discrepancies arose out of inaccurate posi- 
tioning of the probe especially far downstream of the 
nozzles. Many measurements were repeated with dif- 
ferent wires. 


Results and Discussion 


Figs. 3-19 present a summary of the data obtained.j 
Both the mean velocity U’ at a point and the root 
mean square of the axial component of the velocity 
fluctuation «’ are nondimensionalized with respect 
to the mean exit velocity U; so that the turbulence 
traverses are essentially plots of the absolute value of 
the turbulent velocities. 


The Ejector 


The mean velocity profile at the entrance (x/d = 
0.32, Fig. 4) shows that the secondary air has acquired 
a uniform velocity equal to */j) of the jet exit velocity 
and the corresponding turbulence profile shows that no 
mixing has yet taken place. The spike approximately 
located at y = 0.6 in. is probably due in large part to 
the shedding of the turbulent boundary layer from the 
nozzle. For comparison, the mean velocity profile and 
turbulence profiles are given (Fig. 3) for the jet without 
an ejector. At the second station (x/d = 1.68, Fig. 
5) the effect of the ejector on the turbulence levels is 
already apparent. The turbulence within the ejector 
rises to values which are of the order of 30 per cent lower 
than those which prevail in the corresponding sections 
of the free jet. The maximum has also been displaced 
away from the center of the jet. Succeeding traverses 
are similar. Within the ejector the highest turbulence 
level is reached four nozzle diameters downstream of the 
exit. Downstream of that point, turbulence levels de- 
crease smoothly until the end of the ejector is reached. 
In the corresponding free jet the maximum turbulence 
level is reached approximately seven diameters down- 
stream of the exit and it is considerably higher; the 
maximum turbulence within the ejector is approxi- 
mately 68 per cent of what it is in the free jet. Strik- 
ingly, this is approximately the amount by which the 
lateral jump in the mean axial velocity has been reduced 
at the nozzle exit by the presence of the ejector. 

At the ejector exit (Fig. 7) the turbulence level is 
still lower than it is at the same section in the free jet, 
this in spite of the fact that by then the core velocity 
of the free jet has decayed to values which are noticeably 
lower than the centerline velocity in the ejector. 

Downstream of the ejector (Figs. 8, 9, and 10), the 
stream behaves somewhat as a free jet. The centerline 


{ For additional data, see reference 1. 
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velocity decreases little at first, and the turbulence 
levels rise to a value which is roughly the same as that 
for a free jet if one takes into account the reduced 
velocity of the “‘new”’ jet. For the ejector used in this 
experiment, the maximums reached within and down- 
stream of the ejector were very nearly equal. 

To point out the development of turbulence with 
downstream distance, the maximum turbulence occur- 
ring at each cross section has been plotted against 
downstream distance on Fig. 11. 

To summarize the experiment: It has been verified 
that an ejector depresses the turbulence of a round jet. 
It can be argued that, on the other hand, turbulence 
levels have risen downstream of the ejector. This is 
true, and is a consequence of the inability of a jet plus 
ejector to entrain as much air in a given length as a 
simple jet. But this rise of turbulence levels far down- 
stream is likely to be unimportant. It has been pre- 
sented as plausible that the intensity of radiated noise 
was a function of a high power of the turbulence level. 
This means that everything else remaining the same, 
most of the noise comes from regions of highest tur- 
bulence and very little from other regions. Thus, a 
reduction of the maximum turbulence occurring any- 
where in the jet would appear to be the most effective 
way to suppress its noise. To illustrate the point, 
assuming that the only way in which the ejector has 
changed the jet flow is by modifying the turbulence 
levels and using the tentative relation p? ~ w’s, we have 
plotted on Fig. 12 a comparison of the x distribution of 
noise source strength on the basis of Fig. 11. 

In view of the experiment reported above, it seems 
at least plausible that, as was presumed, turbulence 
levels are decreased because the effective velocity of 
the jet with respect to its surroundings is decreased. 
We have seen that both the lateral velocity jump at the 
nozzle exit and the maximum turbulence levels were 
reduced by roughly 30 percent. In order to determine 
whether an ejector does not inhibit turbulence by re- 
straining the lateral motion of the jet, two traverses 
were obtained at x/d = 44 (see Fig. 12). One was 
taken in a free round jet, the other in a jet which was 
trailed by two (cruciform) splitter plates. The latter 
traverse bisected one of the right angles formed by the 
splitter plates. The results show that the turbulence 
levels are not decreased by the splitter plates even 
though the lateral motion of the jet as a whole must be 
inhibited by them. 


The Corrugated Nozzle 


Figs. 14-19 show the results of hot wire traverses at 
eight stations downstream of the corrugated nozzle 
and at two azimuths: @ = 0 traverses bisect the pri- 
mary lobe, @ = 45° traverses bisect the secondary 
lobe, as shown on Fig. 14. 

At the exit of the nozzle, the mean velocity is uniform 
across the nozzle and the turbulence levels are quite 
low except at the sharp boundary of the jet where they 
are again representative of the turbulence in the bound- 
ary layer shed from the nozzle walls. 
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The first mean velocity and turbulent traverses in 
the secondary lobe (Fig. 14, 0 = 45°) reveal that second- 
ary air is entrained before it is mixed in a manner 
analogous to that discussed in connection with the 
ejector. The turbulence levels reach approximately 
the same value for the two azimuths at the first two 
stations. Downstream of the station x/d = 0.8; 
however, the traverses along the primary lobe show 
consistently higher turbulence levels than those along 
the secondary lobes. The maximum turbulence is 
reached for x/d = 2.75 (@ = 45°) and x/d = 3.25 (6 = 
0). Further downstream, the difference between the 
traverses at the two azimuths seem to decrease. Fig. 
19 is a plot of the maximum turbulence at a section asa 
function of downstream distance for the two azimuths. 
It exhibits two features of interest. The first has 
already been noted: wherever secondary air has been 
accelerated before mixing—i.e., in the secondary 
lobes—the turbulence levels are lower. To be sure, 
the reduction of the effective velocity jump and of the 
turbulence levels is smaller than it is for the ejector 
(compare Figs. 19 and 11). On the other hand, the 
corrugated nozzle used was rather shallow in comparison 
with those which have been found most successful in 
silencing jets. In addition, in a deeply corrugated 
nozzle, the fraction of the total perimeter of the jet 
which is affected by the accelerated secondary stream 
is larger. It seems that the lobes of a corrugated noz- 
zle ‘interfere’ with each other in that they compete 
for an insufficient local supply of secondary air in be- 
tween them and thus cause a pressure drop which in- 
duces an axial velocity of the secondary air. In other 
words, secondary air is pumped from upstream in be- 
tween the lobes of the jet. Hence, the effective velocity 
of part of the jet (the part contiguous to the secondary 
lobe) is reduced and this is probably the reason why 
turbulence levels in that region are observed to be lower 
than they would be in a round jet. 

The second feature of the results summarized by Fig. 
19 is that the corrugations have effectively changed the 
length scales of the mixing process. For instance, in a 
free round jet the mean velocity at the core begins to 
decay approximately 6 diameters downstream of the 
nozzle exit and the maximum turbulence is reached 
approximately 7 diameters downstream of the nozzle 
exit. In the corrugated nozzle tested, the centerline 
velocity decay starts 1.6 (original) diameters down- 
stream of the nozzle exit and the maximum turbulence 
is reached approximately 3 diameters downstream. 
For a nozzle with a large number of corrugations, this 
effect is no doubt a great deal more pronounced. The 
success of very short ejectors used in conjunction with 
corrugated nozzles is due to this change of scale. 

In conclusion, the experiments reported above 
provide evidence that a subsonic jet modified so as to 
radiate less aerodynamic noise for a given thrust and 
mass flow exhibits lower maximum turbulence levels 
than an unmodified jet. Experimental data further 
suggest a simple explanation for these results. 


(Continued on page 730) 
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Summary 


Viscous self-induced pressures on 3° semivertex angle cones 
were measured over the Mach Number range 3.7 < M4 < 5.7, and 
for values of the viscous interaction parameter in the range 0.5 < 
X- < 2.3. The data were found to be in good agreement with re- 
sults obtained by Talbot on 5° cones in the range 3.7 < MM, < 4.1, 
0.9 < xX. < 3.6. All these data were correlated reasonably well 
by the viscous interaction parameter. 

A new method for calculating self-induced pressures is presented 
which takes into account the interaction between boundary- 
layer growth and the inviscid flow field at the outer edge of the 
boundary layer. Pressures calculated by this method were only 
10 to 20 per cent higher than the measured values, 


Symbols 


C = Chapman-Rubesin factor in relation (u/u2) = C(T/T2) 
d = diameter of pressure orifice 
K, = similarity parameter, 1/0, 
K» = similarity parameter, 1,6. 


M = Mach Number 

p = static pressure 

ry. = cone radius 

Re = Reynolds Number 

t = bluntness of cone tip or thickness of leading edge of plate 
T = gas temperature 

u = gas velocity 

x = distance measured from vertex along cone surface 
y = distance normal to cone surface 

y = specific heats ratio 

6* = boundary-layer displacement thickness, 


ff, 1 (pu/p2tt2)\dy 


6 = boundary-layer thickness 
u = absolute viscosity 

vy = kinematic viscosity 

p = gas density 

o@ = Prandtl Number 


A = mean free path 

65 = streamline inclination at outer edge of boundary layer 
6. = cone se-nivertex angle 

02 = 0. + % 

X- = viscous interaction parameter, 


Presented at the Aerodynamics—I Session, IAS 27th Annual 
Meeting, N.Y., Jan. 26-29, 1959. 

+ This work was sponsored by the National Advisory Commit- 
tee for Aeronautics. An extended version is contained in NACA 
TN 4327. 

* Associate Professor of Aeronautical Sciences, University of 
California, Berkeley. 

** Professor of Engineering, University of Nagoya, Japan. 
Visiting Professor of Aeronautical Sciences, University of Cali- 
fornia, Berkeley, July-December 1956. 

*** Research Engineer, University of Michigan, Ann Arbor. 


Subscripts 
= free-stream conditions 

2 = conditions at outer edge of cone boundary layer, taken 
to be functions of x 

c¢ = ideal flow conditions along surface of cone obtained 
for x © 

w = wall conditions 

aw = adiabatic wall 


Examples of the notation used for the various Reynolds Numbers 
are: Res/in. = t2/m, Rey, = mt/n, Rez, = ete. 


Introduction 


6 lev FLUID-DYNAMIC and thermodynamic phenomena 
associated with flight at hypersonic speeds have 
been the subject of intensive research in recent years. 
In this research one problem that has received consid- 
erable attention is the problem of the ‘‘self-induced pres- 
sure’ effect, which is one aspect of a broader class of 
phenomena which can be described as ‘“‘viscous inter- 
action’ phenomena. 

Although the fundamental mechanism responsible for 
self-induced pressures is well understood, the analysis 
of the effect is rather complicated. The magnitude 
of the self-induced pressure is directly proportional to 
the rate of growth of the boundary layer. (See Fig. 1.) 
However, the growth of the boundary layer is deter- 
mined by the pressure, Mach Number, etc., in the flow 
at the outer edge of the layer, and the values of these 
quantities depend on the magnitude of the displace- 
ment effect. It is seen, therefore, that we have to deal 
with a complex interaction phenomenon in which the 
boundary-layer “‘history’’ plays an important role. 
It will also be recognized that the phenomenon is, 
with regard to magnitude, more significant for thin 
bodies such as flat plates and slender cones than for 
thick bodies, since for thin bodies the changes in ef- 
fective geometry due to boundary-layer growth will be 
proportionately larger. 
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Fic. 1. Schematic of viscous flow over a cone. 
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DIMENSIONS INCHES 
[a |-9.50/0s00] 0.75 |1.25 | 2.00] 2.75 | 3.50] 4.00 | 4.50] 
rte sooecs [1a 75 [0.730] 0.75 | 3.00 | 4.50 | 6.00] 6.50 | 7.00] 
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MATL: ALL BRASS EXCEPT FOR ', 0.0. SOFT COPPER TUBE 


Fic. 2. Model specifications. 


The particular problem considered in this paper is 
the self-induced pressure effect on a slender cone. In 
the first part of the paper new data are presented for 
pressures on 3° semivertex angle cones in the range of 
flow conditions 3.7 < M, < 5.8; 650 < Re,, < 32,000. 
Theories used for comparison with experiment are 
presented in the Appendixes. 


Previous Experimental Results 


The only data which have been available up to now 
for the viscous interacting air flow over a cone are the 
experiments of Talbot! and Baldwin.? Talbot's results 
were obtained in the Low Density Wind Tunnel of the 
University of California, the same wind tunnel used 
for the present work. His tests were carried out on 5° 
semivertex angle cones, in the Mach Number range 
3.7 < M, < 4.1, and at Reynolds Numbers which cor- 
responded to the viscous interaction parameter range 
0.91 < X. < 3.54. Baldwin’s data were obtained in 
the GALCIT 5” X& 5” Hypersonic Wind Tunnel, Leg 
No. 1. His data were also for a 5° cone, at MJ, = 
5.8, and over the viscous interaction parameter range 
EG: 


Description of Present Experiments 
Wind Tunnel 


The experiments were conducted in the No. 4 Low 
Density Wind Tunnel of the Low Pressures Project of 
the University of California. This wind tunnel, which 


TABLE 1 
Flow Conditions of Tests, 7; = 300°K. 

pr | Pr 
M, (uHeg) Re,/in. | (0. = 3.03°) 
7 49.8 880 | 1.106 
91 73.0 1,510 | 1.114 
.97 85.1 1,860 eg 
108.8 21540 1.120 
47 66.5 4,530 1.207 

113.3 8,980 1.223 
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is an open-jet continuous-flow type employing axially 
symmetric nozzles, is described in detail elsewhere.’ 
Two nozzles were used in the tests: the No. 8 nominal 
Mach 4 nozzle which produces flows in the range 3.7 < 
M, < 4.1, 90 < Re,/in. < 3,600, and the No. 9 nominal 
Mach 6 nozzle, which produces flows in the range 5.5 < 
M, < 5.8, 4,000 < Re,/in. < 9,000. Actual values of 
the flow parameters obtained in the test are listed in 
Table 1. 


Models 


All the cones tested were of 3° semivertex angle. 
Two sets of models were used, each set consisting of 
seven cones (see Fig. 2). The Type A cones had base 
diameters of 0.500 in., and the Type B cones had base 
diameters of 0.750 in. The longer Type B models were 
designed primarily to investigate the influence on the 
cone surface pressure of the expansion wave generated 
at the juncture of the conical surface and the cylindrical 
afterbody. Alignment of the cones in the flow was 
accomplished by adjustable set screws in the base 
support. 

Subsequent to its use in the pressure measurement 
tests, Model B-7 was fitted with four copper-constantan 
thermocouples soldered in the cone surface at 2, 3, 4, 
and 5 in. from the vertex. This model provided in- 
formation on the wall temperature of the cones, which 
was required for the boundary-layer calculations. 

Some difficulty was encountered in producing models 
with sharp tins. The method of fabrication which was 
found most satisfactory was an acid etching process. 
After the cones had been machined to nearly their 
final dimensions, the tip regions were etched by a flat 
tool covered with a thin film of nitric acid. This 
method produced tips with diameters less than 0.001 
in. However, the etching was not completely uniform. 
In the etched regions of the models, which extend back 
from the vertices about 0.4 in., some local variations in 
cone angle of several degrees were observed. Further 
back from the vertices all cone angles were found to be 
3.08° + 0.08°. 
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Fic. 3(a). 3° cone data, M ~ 4, 
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Instrumentation 


Cone surface pressures were measured with a tem- 
perature-regulated thermistor manometer. The least 
count (0.1 mv) of the potentiometer used to measure 
the bridge unbalance of the thermistor measuring cir- 
cuit corresponded to a pressure increment of about 
0.02 micron Hg. Before each test the thermistor was 
calibrated statically against a precision McLeod gage. 
Analysis of the thermistor calibration data yielded a 
probable error in absolute pressure of about 1 per cent. 

Wind-tunnel stagnation pressures were measured 
with a mercury manometer. Impact pressures were 
measured with a butyl phthalate oil manometer. Both 
manometers were equipped with magnifying optics 
which made it possible to locate the menisci to within 
0.001 in. 


Nozzle Calibration 


Mach Numbers and static pressures in the test re- 
gions of the nozzles were determined by measuring 
stagnation and impact pressures, assuming the flow 
to be isentropic. For each flow condition of the tests, 
an axial impact pressure survey was made to determine 
the Mach Number and static pressure variations in the 
regions occupied by the models. 

The tunnel-empty traverses in the No. 8 nozzle 
revealed a region about 8 in. in axial extent over which 
the Mach Number variation was less than 2 per cent 
and the static pressure variation less than 10 per cent. 
In the No. 9 nozzle the axial extent of the region 
over which the Mach Number and static pressure vari- 
ations were less than these values was about 3.5 in. 


Procedure 


In the tests reported, all models were positioned so 
that their vertices were located at the same point in 
the flow. Correction was made for the axial gradients. 
Other tests were made with models positioned so that 
their respective pressure orifices were at the same axial 
location in the stream. The surface pressures obtained 
by this second method agreed quite well with those ob- 
tained by the first method. The correction for axial 
gradients in the stream was accomplished simply by 
using the local Mach Number (from tunnel-empty 
measurements) in the determination of the inviscid 
cone pressure p,. It was found that this procedure gave 
consistent results provided the static pressure in the 
flow at the point where orifices were located did not 
differ by more than about 10 per cent from that at the 
vertex of the model. Consequently, in the No. 9 
nozzle where the length of usable flow was about 3.5 in., 
models for which x = 3.5 in. were not used in the final 


tests. No difficulties arose on this score in the No. 8 
nozzle. However, it was observed that the expansion 


wave reflection of the bow shock (the reflection occur- 
ring in the region of strong density gradient where the 
isentropic core merged with the nozzle boundary layer) 
affected the pressures measured on models B-6 and B-7. 
For this reason data are not reported at / ~ 4 for these 
models. 
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Fic. 3(b). 3° cone data, M ~ 4. 
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Fic. 3(d). 3° cone data. M ~ 4. 
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Fic. 4(a). 3° cone data, M ~ 6. 
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Fic. 4(b). 3° cone data, MW ~ 6. 


The temperature measurements made with Model B-7 
showed the surface of the cone to be isothermal. The 
apparent recovery factor (based on inviscid flow condi- 
tions I/,, T, behind the conical shock) was about 0.89; 
the increase over the theoretical value of 0.85 for the 
laminar boundary layer is due to heat conduction 
through the support from the model mounting, which 
was at essentially stagnation temperature. 


Results 


Results of the tests are given in terms of the induced 
pressure increment p2./p, — 1, where fp» is the measured 
cone pressure and /p, the inviscid Taylor-Maccoll value. 
The inviscid flow values for the 3° cones were calcu- 
lated from Van Dyke’s second order theory,‘ since 
interpolation in the Kopal Tables®> between 0° and 5° 
was not sufficiently accurate. The data at different 
flow conditions are shown plotted versus x in Figs. 3 and 
4, and in Fig. 6 versus the hypersonic similarity param- 
eter x. = M,5(C/Re,,)"? where M, and Re,, are the 
Mach Number and Reynolds Number based on ideal 
Taylor-Maccoll flow conditions. 

Numerical values for the induced pressure increment 
are in the range 0.06—0.30. The actual measured cone 
pressures varied between about 60 and 170 microns Hg. 
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An error of 1 per cent in the measured pressure is 
equivalent to an error in the induced pressure increment 
of 5 per cent or more for most of the data, assuming the 
values of p, to be exact, and from this it is estimated 
that the overall probable errors in p2/p, — 1 are between 
5 and 15 per cent. Free-stream Mach Numbers are 
accurate to about 1 per cent; free-stream Reynolds 
Numbers are accurate to about 5 per cent. 

In addition to the results of the present tests, two 
sets of data taken from reference 1 are shown in Fig. 5. 
In Fig. 6 all of the data from reference 1 on models A-1 
through A-7 are plotted, and also Baldwin’s 5° cone 
data are represented by a single line. 


Discussion 


Experimental Results 

One conclusion which can be drawn from an exami- 
nation of Figs. 3 and 4 is that the effect of the shoulder 
expansion did not extend far enough upstream in the 
cone boundary layer to influence the cone surface pres- 
sures, since the data obtained with the B models agree, 
within experimental scatter, with those obtained with 
the A models. We were not able to determine the ex- 
tent of the region on the cone which is influenced by the 
shoulder expansion, because the reflection of the bow 
shock wave back onto the models obscured the effect. 
Models B-6 and B-7, which were designed to measure 
upstream influence, were those most affected by the 
reflected wave. In Fig. 3(a) it will be noted that the 
pressure orifices on model B-5 were also within the zone 
influenced by the reflection. 

The scatter in the data is probably due to a combi- 
nation of experimental error, imperfections in the orifices, 
and inaccuracies in the cone angle in the tip regions of 
the models. One may note that the reproducibility of 
the data was quite good, as evidenced by the com- 
parisons shown in Figs. 4(a) and 4(b) between different 
sets of measurements made with both the A and the B 
models. 

It can be seen from Fig. 6 that the parameter x, pro- 
vides a fairly good correlation for all of the data ob- 
tained in the University of California Low Density 
Wind Tunnel. The Mach 4 data for the 3° and 5° 
cones agree quite well; the Mach 6 data are slightly 
lower. 

it is also seen from Fig. 6 that the induced pressure 
increments found by Baldwin are higher by a factor of 
about 2 than those obtained here. One suggestion 
which has been advanced is that the differences may be 
due in part to the influence of tip bluntness. It is true 
that the Reynolds Numbers based on tip diameter 
were higher in Baldwin’s experiments than in ours; 
Baldwin’s were mostly in the range 65 < Re, < 230, 
whereas all the Low Density Wind Tunnel data corre- 
spond to Re;, <9. However, it seems unlikely that tip 


bluntness could account for much of the difference. 
The experiments on flat plates* * indicate that below 
about Re,, = 80 the plate can be considered as “‘sharp,” 
and the effect of tip bluntness has been shown to be 
much less pronounced for cones than for flat plates.® 
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Comparisons Between Theory and ‘Experiment 


Three methods for calculating self-induced pressure 
have been employed in this report. The first, a new 
method devised by us, is presented and discussed 
in Appendix (A). The second, a method proposed by 
Probstein, is reviewed in Appendix (B), as is the third, 
which we have called the ‘“7C.."’ method. 

Figs. 3, 4, and 5 show that the present method for 
calculating self-induced pressures generally overesti- 
mates the data obtained in the Low Density Wind 
Tunnel by about 10-20 per cent. In contrast, the 
TC. method and the Probstein second-order theory 
both give values greater by about a factor of 2 than the 
experimental results. The better agreement obtained 
with the present method is not surprising. Of the 
three methods it is the only one which accounts in even 
an approximate way for the true interaction effect, 
wherein the changes in the external flow due to the 
presence of the boundary layer feed back into the layer 
and alter its rate of growth. It seems likely that the 
discrepancies remaining between experiment and theory 
may be due mainly to the transverse curvature effect. 
[See Appendix (A). ] 

In Fig. 6 the results of the present theory are repre- 
sented by a single straight line. (Actually, the indi- 
vidual curves of Figs. 3, 4, and 5 when plotted against 
x- deviated about +2-4 per cent from this mean line. 
The deviations seemed not to follow any particular 
trend, and were probably due mainly to accuracies 
introduced in the graphical parts of the analysis.) 
Again, in Fig. 6, the good agreement between the pres- 
ent theory and the Low Density Wind Tunnel data is 
clearly evident. Induced pressures have also been 
measured on 5° cones in helium” in the range JJ; = 16 
to 18, and 0.8 < xX, < 1.6. A best fit through the higher 
points of these data is (p2 — p,)/pe ~ 1 + 0.25 X,, and 
if we correct to a y = 1.4 gas by the theoretical weak- 
interaction factor!’ y(y — 1), we obtain (po — p.)/p. * 
1 + 0.13 X., which may be compared with the best fit 
for the present data, (p2 — p,)/pe ~ 1 + 0.12 Xe. 


Hole Size Effect 

It has been shown by Talbot! and by Rayle,® and 
others, that the apparent pressure sensed by a static 
pressure orifice increases with the diameter of the orifice. 
The phenomenon is due to mixing between the stream 
passing over the surface and the fluid confined within 
the orifice and pressure tubulation; the momentum 
transferred by the mixing sets up currents in the fluid 
within the orifice which give rise to the increase in 
pressure. 

Ideally, a static pressure orifice should be as small as 
possible, both to minimize this hole size effect and to 
provide a truly localized pressure measurement. How- 
ever, in rarefied gas flow if a pressure orifice is made 
small enough one encounters another effect, known as 
thermal transpiration, which can also result in errors 
in pressure measurement. Thermal transpiration oc- 
curs, for example, when an orifice whose diameter is 
small compared to the mean free path separates two 


regions of gas at different temperatures.'° In this 
case, the pressure ratio is given by 


= (1) 


The static pressure orifices used in the present experi- 
ments were 0.010 in. in diameter. For an orifice of this 
size, the pressure increment due to momentum mixing 
is completely negligible. However, there is the possi- 
bility that thermal transpiration effects may be impor- 
tant, since the boundary layer is a region of strong tem- 
perature gradient, and many of the molecules which 
enter the orifice from the gas stream come from regions 
in the boundary layer which are at temperatures differ- 
ent from the gas within the orifice. 

We can make a rough estimate of the magnitude of 
the thermal transpiration effect in the following way. 
We assume the boundary-layer characteristics to be 
given with sufficient accuracy by Howarth’s analysis.!! 
For an insulated cone, with o = 1 and w/w. = T/T», 


= (5.0/V3)(1 + 0.08M22) V vox /us (2) 


T/Tow = {1 + 0.2M2[1 — (u?/u?)]}/(1 + 0.2M2?) 
(3) 
Now, let us also assume that the velocity distribution 


in the boundary layer is linear in y. Then the tem- 
perature 7 at a distance A, from the wall is given by 


T/T = 
{1 + 0.2Me2[1 — (A,2/62)]}/(1 + (4) 


0 = 5°, mM, = 3.70, Re, = 930/in, 


PRESENT THEORY 
PROBSTEIN - 15' ORDER 
PROBSTEIN - 2°° ORDER 
3 TC@ METHOD 
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Fic. 5. 5° cone data, M ~ 4. 
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Fic. 6. Induced pressure increment vs. viscous interaction 
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We identify A, with the mean free path of the gas at 


the wall. 
As a specific example, let us take the following condi- 


tions: 


M, = 5.5 x = lin. 
Re, = 10,000/in. Tw = 270°K. 
pe = 170 pHg 

For these conditions A, = 0.011 in., 6 ~ 0.099 in. 


Then from Eq. (1) we find p,/py. = 0.994, where we 
identify p, as the pressure which we are attempting to 
measure and p, as the pressure within the orifice which 
is presumably in error because of the thermal tran- 
spiration effect. It is seen that for this particular case 
the error is about 1/2 per cent. For the worst condi- 
tions, the error is found to be about 2 per cent. Actu- 
ally, this analysis greatly overestimates the effect, since 
Eq. (1) is true only for A,/d > 1. For the experiments 
A,/d ~ 1, and in this range the pressure increment is 
less than 10 per cent of what it is in free molecule flow. 
One may conclude, therefore, that the pressures meas- 
ured in the experiments were true static pressures, essen- 
tially uninfluenced by either momentum mixing or 
thermal transpiration hole size effects. 


Appendix (A) 
“‘Modified Tangent-Cone’’ Method 


The assumptions usually made in the tangent-cone 
analysis for self-induced pressures are: (a) a region 
of inviscid flow exists between the outer edge of the 
boundary layer and the shock wave, and (b) the flow 
parameters, such as Mach Number, pressure, etc., at 
the edge of the boundary layer (subscript 2 in Fig. 1) 
can be obtained to satisfactory approximation by the 
tangent-cone (7C) method. The 7C method consists 
of relating the local ow parameters on a body to the 
undisturbed flow ahead of the shock wave through 
conical flow theory (e.g., Taylor-Maccoll values), but 


6*/x = VC/Re,( (4/2) — 
tan 6; = (dé*/dx) = (1/2)-(6*/x) 
where = 1 + — 1)M22/2 (A-2) 
and the Chapman-Rubesin C is defined by 
C = (A-3) 


The viscosity u’ is evaluated at the intermediate 7’ 
given by 


T'/T. = (Tx/T.) — 0.4680" [(T~ — Tax)/Te] 


0.2730[(y — (A-4) 


However, instead of using the Sutherland law, as recom- 
mended by Monaghan, the Bromley-Wilke'® values 
are used for u’/u,, since the Sutherland law is less ac- 
curate at the low free-stream temperatures encountered 
in our tests. The actual wall temperature 7, is deter- 
mined by experiment. 
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using the local body inclination as the effective cone 
angle. In the case of self-induced pressures on a cone 
as shown in Fig. 1, the local effective cone angle is 
taken to be 62, the sum of the cone angle @, and the angle 
6; = tan~'(dé*/dx). 

Two computations are involved in utilizing the 7C 
method for evaluation of self-induced pressures. _ First, 
the inviscid flow values must be obtained, for given 
effective cone angle 62, either from the exact solution 
(Kopal’'s Tables) or by one of several approximate 
methods!® '* which are available. Second, the bound- 
ary-layer displacement thickness must be evaluated, as 
a function of position along the cone surface, and here 
again several methods of varying accuracy and com- 
plexity are available. It will be noticed that the two 
computations are not independent. We are dealing 
with an interaction phenomenon—the boundary-layer 
growth determines the inviscid flow values at the outer 
edge of the boundary layer, but at the same time the 
rate of growth of the boundary layer is determined by 
these inviscid flow values. An accurate application 
of the 7C method must include this interaction effect. 

For the computations of the inviscid flow values, 
the exact Taylor-Maccoll results as computed by Kopal 
areemployed. There are several supersonic and hyper- 
sonic flow approximations in analytic form which are 
accurate over different ranges of the similarity param- 
eter K, = 1,62, but no single one is sufficiently accurate 
over the entire range of AK» encountered in the tests 
reported here. The boundary-layer displacement 
thickness is calculated from the approximate formula 
of Monaghan,'* which includes the effect of Prandtl 
Number and isothermal-wall heat transfer. It does 
not include the effects of pressure gradient or trans- 
verse curvature. To account in some measure for the 
variation in the external flow qualities along the outer 
edge of the boundary layer, local values of the external 
flow parameters are used to calculate \/. and Ke,,. 

Monaghan’s result, to which the Mangler’ trans- 
formation correction has been applied, is 


— {1 + — 
\ 


(A-1) 


For a given free-stream condition (/;, fp, etc.), a 
set of values of 62 are chosen, and for each of these values 
a range of values of @; are computed, using local free- 
stream conditions determined by the tangent-cone 
method and Kopal’s Tables. The angle 6; is a function 
of x, for each selected value of 62. Since the value of x 
appropriate to a particular value of 6. is that for which 
6. — 0; = 6., by cross plotting one can obtain 6) as a 
function of x, and, hence, p2/p, as a function of x. 


Accuracy of the Method 


The accuracy of the TC calculation of the inviscid 
flow has been examined by Ehret!’ and Lees,'* by com- 
paring pressure distributions on pointed ogives with 
exact values from the method of characteristics. 
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HYPERSONIC VISCOUS FLOW OVER SLENDER CONES 


TC method yields surface pressures which are slightly 
higher than the exact values, the difference depending 
on the distance from the vertex of the ogive. At the 
vertex the two methods, of course, give identical re- 
sults; farther back the deviation may be of the order of a 
few per cent. It is found that the 7C method also 
overestimates the pressure for blunt power-law bodies. '* 

In the calculation for self-induced pressures the effect 
of boundary-layer growth on the external inviscid flow 
is approximated by increasing the effective local cone 
angle by @;. Order of magnitude arguments concerning 
the accuracy of this approximation have been given by 
Lees and Probstein.*”” The error involved in replacing 
the actual streamline inclination in the external flow 
by 4; is estimated to be of order (6/x)*, where 6 is the 
boundary-layer thickness. It also turns out that the 
neglect of pressure gradient across the boundary layer 
is also justified provided (6/x)? is small. 

The variation in external flow properties along the 
outer edge of the boundary layer is partially taken into 
account by using local values for Mo, Reo, po, ete. But 
it will be noted that the expression for dé*/dx [Eq. 
(A-1)] is only approximate, since the terms involving 
dM./dx, dT2/dx, and dRe./dx have been neglected. It 
turns out that for the present calculations these terms 
contribute an increment of about 5 per cent at most 
to #;, and their neglect is not serious. 

Two important effects which have not been included 
in the boundary-layer analysis are the effect of trans- 
verse curvature and the direct effect of the self-induced 
pressure gradient on the density and velocity distribu- 
tions within the boundary layer. Both the transverse 
curvature and the pressure gradient tend to thin the 
boundary layer, and thus result in smaller values for 
the induced pressure increment. The _ transverse 
curvature effect has been studied by Probstein and 
Elliott.*!_ Probstein®?? concludes from his analysis, 
which is valid for small 6*/r,, that transverse curvature 
does not appreciably alter the boundary-layer displace- 
ment thickness. However, much of our data were in 
the range 6*/r, = 1 to 3, and for these values the trans- 
verse curvature effect almost certainly cannot be neg- 
lected. An estimate of transverse curvature effect 
in very simple form is given by Hill, et al.,** in the 
form 

(dé*/dx),. = (dé*/dx) Vi + 25* re 

where 6* is the displacement thickness calculated by 
theory, neglecting transverse curvature, and 6,,* the 
actual displacement thickness. This correction is 
probably not valid for values of 6*/r, as large as those 
of the present tests, but we may note that for 6*/r, = 
1, a reduction of about 40 per cent in the induced pres- 
sure is predicted. 


Appendix (B) 
Probstein’s Analysis for the Self-Induced Pressure Fffect 


Probstein®’ considers a Taylor series expansion of the 
surface pressure in the form 


(pz — Pe)/Pe = (br/Pe) [O(Pe/ pr) * + 
(1/2!) [0°(p./p1) O52 +... 


Now, for K, < 1, which is the range encountered in the 
present tests, the Lees hypersonic approximation is not 
accurate, so that the values for the derivatives in (B-1) 
obtained by Probstein cannot be used. Probstein 
suggests that for A, < 1 the ratio p,/p, be evaluated 
from the Kopal Tables, but that the derivatives be 
evaluated from the von Karman** slender body result. 
However, this result is also not sufficiently accurate. 
Guided by Van Dyke’s suggestions® for combined 
supersonic-hypersonic similarity, we found that a for- 
mula of the form 


(B-1) 


Pr = 1 + In M2 — 1) (B-2) 


could be made to fit the Kopal values and the Van Dyke 
second-order theory values for p,/p; over the range 
0.14 < — 1 < 0.3; 0 < & < 0.13 radians, 
The constants found were A; = 1.52, Ay = 2.85. Ex- 
pression (B-2) was used to calculate the derivatives in 
Eq. (B-1), and expression (A-1) was used to evaluate 
6; and @;", except that inviscid flow values /,, Re,, and 
T., were used instead of the local values J/>, Res, and 

We found that it was necessary to include at least 
two terms in the series. For example, with 1, = 


3.70, 0. = 5°, 


(pi Pc) [O(p, pe), 06 =0, = 3.5 


B-3 
Pc) = 17 (B-3) 


We also performed a numerical differentiation of a 
curve constructed from cross-plots of the Kopal entries 
and the 3° values calculated by second-order theory, 
and obtained the values 3.8 and 38 for the above deriv- 
atives. The first of these values may be more accurate 
than that given in Eq. (B-3). Not much accuracy can 
be claimed for either set, however. 

The results of the Probstein analysis are shown in 
Figs. 4(b) and 5. The curves entitled “Ist Order” 
were obtained using only the first term in the series 
(B-1), those entitled ‘‘2nd Order’’ using two terms. It 
can be seen that for small angle cones (i.e., small A,) 
the convergence of the series is slow. This slow con- 
vergence is already evident at K, = 1, as can be seen 
from examination of the functions presented in refer- 
ence 22. 

As a check on the Probstein analysis, we also evalu- 
ated two induced pressure distributions by what we 
have called the “7C.." Method. In this method 
the boundary-layer slope 6; is calculated from the in- 
viscid flow values 7., M,, and Re,, rather than the 
local values of these quantities, but the tangent cone 
method with the Kopal values rather than the Prob- 
stein series is used to evaluate the pressure. The 
TC.. method, of course, gives the values that the Prob- 
stein method should converge to. As is evident in Figs. 
4(b) and 5, the use of two terms in the Probstein series 
provides a fairly good approximation. But Probstein’s 
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method is not very useful for K, < 1, because it is diffi- 
cult to obtain accurate values for the required deriv- 
atives. 
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Some Effects of Sound-Reduction Devices on a Turbulent Jet 


(Continued from page 722) 


One cannot, on the strength of the present experi- 
ment alone, determine whether the turbulence values 
u’/U; yield a sufficiently complete working index of 
the noise radiating efficiency of a jet. Still, reasons 
have been adduced to hope that the evaluation of the 
noise generating strength of a jet by its turbulence 
levels will prove adequate for engineering purposes. 
To establish the point, the type of measurements 
reported here should be made on a number of diiferent 
nozzles with known noise generating efficiency. 

What has been said in no way obviates the need for 
a more informative (and more difficult) study which 
might relate a more exact measure of the noise source 
strength with some integral property of the turbulent 
region 
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Heat Transfer in Separated Flows 


HOWARD K. LARSON* 


Ames Research Center, 


Summary 


Results of an experimental heat-transfer investigation in re- 
gions of separated flow are presented and compared with the 
theoretical analysis of NACA TN 3792. The average heat 
transfer for both laminar and turbulent separated boundary 
layers was found to be from 35 to 50 per cent less than that for 
equivalent attached boundary layers. The overall scope of the 
measurements included Mach Numbers from 0.3 to 4.0 and 
Reynolds Numbers from 165 to4 & 10°. The results for laminar 
boundary layers agree well with the analysis of TN 3792. The 
results for turbulent boundary layers, however, disagree con- 
siderably. Results of velocity and temperature surveys in the 
separated turbulent boundary layer are presented and partially 
explain the discrepancy between the experiments and analysis. 
The maximum local heat-transfer rates were found to occur in 
the reattachment region of the separated boundary layers investi- 
gated. The effect of transition on heat transfer in the separated 
laminar boundary layers is described and data showing effects 
of Mach Number and wall temperature on the transition Reyn- 
olds Number of separated laminar flows are also included. 


Symbols 

Ap, = area of heated boundary layer: 7D¢ for axially sym- 
metric models; 2bf for two-dimensional models, 
sq.ft. 

b = span of two-dimensional models, ft. 

¢, = specific heat of air, B.t.u./Ib. °R. 

D = diameter of axially symmetric models, ft. 

h = average heat-transfer coefficient, h = d(Qtor/Agz,)/dT, 

t = length of heated boundary layer, ft. 

{’ = length of separated boundary layer, ft. (€’ = ¢ for 
axially symmetric model, {’ == 2/3 ¢ for two-dimen- 
sional model ) 

L = model length, ft. (defined in Fig. 2) 

M = = Mach Number 

Qtot = total heat-transfer rate to heated portion of model, 
B.t.u./sec. 

q = local surface heat-transfer rate, B.t.u./sec. ft.* 

q = average surface heat-transfer rate, B.t.u./sec. ft.* 

R, = Reynolds Number, = p,u,L/n, 

R,, = transition Reynolds Number, R,, = potlofl’/po 

S = separation location 

St = average Stanton Number, Si =h/p,u,c, 

T = temperature, °R. 

u = velocity in x direction, ft./sec. 

x = longitudinal distance downstream of separation, ft. 

y = radial distance from model surface, ft. 

= depth of separated region, ft. 
p = density, lb./cu. ft. 
m = coefficient of viscosity, lb./ft. sec. 
Subscripts 

A = attached boundary layer 

aw = adiabatic wall 

e = outer edge of separated or attached boundary layer 


Presented at the Aerodynamics—II Session, IAS 27th Annual 
Meeting, N.Y., Jan. 26-29, 1959. 
* Aeronautical Research Scientist. 
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S = separated boundary layer 
w = wall 
oo = free stream 


Introduction 


A KNOWLEDGE of the effects of flow separation on 
heat transfer is an important consideration in the 
design of a variety of aerodynamic vehicles. In refer- 
ence 1, Dr. D. R. Chapman of Ames Research Center 
presented an analysis of heat transfer in regions of sep- 
arated flow. The results of that analysis indicate that 
separation of a laminar boundary layer reduces the aver- 
age heat transfer whereas separation of the turbulent 
boundary layer can either increase or decrease the 
average heat transfer, depending on the Mach Number. 
Experimental heat-transfer data, appropriate for com- 
parison with the results of the analysis, have not been 
available in the literature. 

In order to provide experimental heat-transfer data 
for comparison with the analytical results, a wind-tunnel 
investigation was initiated. These experiments pro- 
vided heat-transfer data for equivalent attached and 
separated boundary layers, thereby enabling a direct 
comparison of separated-flow heat transfer with at- 
tached-flow heat transfer. Additional observations of 
this investigation also provided an insight into some of 
the details of heat transfer in separated boundary layers. 

This paper presents experimentally determined ratios 
of average heat transfer for separated boundary layers 
to average heat transfer for equivalent attached bound- 
ary layers and compares the results with the calcula- 
Data were obtained for completely 
Since 


tions of reference 1. 
laminar boundary layers in supersonic flow. 
reference | indicated a large effect of Mach Number on 
heat transfer for separated turbulent flows, data for 
turbulent boundary layers were obtained over a range 
of subsonic and supersonic Mach Numbers. Qualita- 
tive measurements of the distribution of local heat- 
transfer rates were obtained and are presented for 
both separated laminar and turbulent flows. Utili- 
zation of fully laminar separated boundary layers 
(transition downstream of reattachment) depends on 
the transition Reynolds Number characteristics of 
these flows. Since aerodynamic heating problems are 
usually encountered at wall temperatures considerably 
lower than adiabatic, some recent measurements of the 
effect of wall temperature on separated boundary- 
layer transition are also included. 

The writer recently became aware of the results pre- 
sented in references 2-6. A short section is added at 
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REAT TACHMENT 


EDGES OF MIXING LAYER ZONE 


Theoretical model. 


Fic. 1. 


the end of this paper for the purpose of comparing the 
present data with the results of those investigations. 


Description of Analysis 


The physical model used in the analysis of reference 
lis presented in Fig. 1. The theory analyzes the thin, 
separated boundary layer with the following assump- 
tions: (a) the thickness of the boundary layer is zero 
at separation, or small compared with the depth of the 
separated region; (b) the length of the reattachment 
zone is small compared with the length of the separated 
boundary layer; (c) the pressure is constant along the 
separated layer; (d) the wall is isothermal; and (e) the 
low-speed air in the reverse-flow region is at wall tem- 
perature when drawn into the separated boundary 
layer at the zero velocity stream line. 

The theory then calculates the average heat transfer 
for the separated boundary layer and compares it with 
the average heat transfer for an attached boundary 
layer with the same constant Mach Number, temper- 
ature, and pressure at the outer edge of the boundary 
layer, and the same constant wall temperature. Com- 
parisons of the experimental results with the calculations 
of the analysis will be presented separately for laminar 
and turbulent flows after first describing the experimen- 
tal method. 


Experimental Method 


Two types of electrically heated models were used 
in the experimental investigation. As shown at the 
left of Fig. 2, the axially symmetric models provided 
cylindrical attached and separated boundary layers of 
approximately the same mean, outer-edge flow prop- 
erties and boundary-layer area. The two-dimensional 
wedge models, at the right, provided two-dimensional 
attached and separated boundary layers with approxi- 
mately the same outer-edge flow properties and area. 
During the tests of both the axially symmetric and two- 
dimensional models, the attached- and separated-flow 
models were tested under nearly identical wind-tunnel 
flow conditions, thereby providing a direct comparison 
of heat transfer for separated boundary layers with heat 
transfer for attached boundary layers. Each of the 
two heated centerbodies of the axially symmetric models 
had eight individually heated segments. The forebody 


and afterbody were not heated and were thermally 
insulated from the heated centerbodies. During the 
tests, power was distributed among the eight segments 
to provide a constant wall temperature. The two- 
dimensional models were constructed of solid copper, 
and power to fore and aft electrical heaters was ad- 
justed to provide a constant wall temperature. In 
order to provide true, two-dimensional data, guard 
heaters of identical cross section were mounted out- 
board on both sides of the models and capped with 
thin insulated plates to prevent crossflow in the sepa- 
rated region. A rear guard heater and the above 
described end guard heaters were removed when 
shadowgraphs of the flow pattern were taken. 

The experimental data were obtained in the follow- 
ing manner: with the power distributed to provide 
constant wall temperature, the total power to the 
model was varied and plotted as a function of wall 
temperature. Fig. 3 presents some typical data plots 
from which the average heat-transfer coefficients, h, 
were calculated. For the purposes of the present tests, 
the average heat-transfer coefficient is defined as ; = 
d(Qtot/Aar)/d7T and is equivalent to the customary 
definition, = — Tow), provided h is 
independent of 7,,. In the course of the present tests, 
h was always found to be independent of 7,, within 
the limits of experimental accuracy. The use of the 
present definition of / does not require a knowledge of 
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the adiabatic wall temperature. During most of the 
runs, power-off data were obtained and the measured 
adiabatic wall temperatures agreed well with extra- 
polations of the Q,o./A xz versus 7), plots to zero heat 
transfer. 

For the purposes of comparing the experimental re- 
sults with the analysis of reference 1, the data will be 
presented as the ratio of the average heat-transfer 
coefficients for separated and attached boundary 
layers, hs/h4, obtained with nearly identical outer- 
edge flow conditions. The Reynolds Number and the 
Stanton Number were calculated using mean flow prop- 
erties at the outer edges of the boundary layers. 


Comparison of Experimental Results for 
Laminar Boundary Layers With Theory 


Axially Symmetric Models 


Data obtained for the axially symmetric models are 
presented in Fig. 4, and compared with the results of 
reference 1. The measurements of fs/h, at M, = 
3.1 and 4.0 agree well with the theory and appear to be 
reasonably independent of Mach Number and Reyn- 
olds Number up to the Reynolds Number at which 
transition begins in the reattachment region. At 
higher Reynolds Numbers, the increased values of 
hs/h, are for transitional separated-flow data com- 
pared with laminar attached-flow data. Separation is 
observed to induce premature transition and at suffi- 
ciently high Reynolds Numbers, R,; > 2 X 10%, an in- 
crease rather than a reduction in average heat transfer 
was measured. 

The interpretation of the agreement between the ex- 
perimental results for completely laminar flow and 
reference 1 was considered open to question since the 
heated area of the separated-flow model was 35 per cent 
less than the heated area of the attached-flow model. 
In order to determine to what extent the reduction in 
heat transfer is due to the reduction in surface area, 
additional experiments were conducted with the two- 
dimensional models described in Fig. 2. The wetted- 
surface area lying beneath the separated boundary layer 
is, for this case, about 10 per cent greater than that for 
the attached boundary layer. 


Two-Dimensional Models 


Data obtained in the study of the two-dimensional 
boundary layers are presented in Fig. 5 and compared 
with theory. These data are for M, = 3.0 and were 
obtained at a free-stream Mach Number of about 4. 
Only the rearward two-thirds chord of the separated- 
flow model was actually exposed to separated flow, and 
the results of reference 1 can be applied only to that 
portion of the surface for comparison of experiment 
With theory. Hence, the full theoretical reduction of 
44 per cent would not be expected. Heat transfer for 
laminar attached boundary layers varies as x~'/*, and 
the theoretical fractional reduction in heat transfer for 


this case would be (4 — hs)/hg = 0.44(VW1 — 
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analysis corresponds therefore to a value of hs/h, = 
0.814. The data again agree well with the theory up 
to the Reynolds Number at which transition begins in 
the reattachment region. The present results then 
indicate that the average heat transfer for separated 
laminar boundary layers is reduced in accordance with 
the calculations of reference 1 and is apparently inde- 
pendent of the heated surface area lying beneath the 
separated boundary layer. 


Comparison of Experimental Results for 
Turbulent Boundary Layers With Theory 


In order to obtain turbulent-boundary-layer data 
for comparison with the analysis, boundary-layer trips 
were fitted to the insulated forebody of the axially 
symmetric models described earlier. Two shadow- 
graphs typical of the turbulent flows investigated are 
presented in Fig. 6. It can be observed in both 
shadowgraphs that transition occurred just upstream 
of separation. This location of transition provided a 
minimum-thickness, turbulent boundary layer at sep- 
aration in accordance with the assumptions of reference 
1. Boundary-layer trip requirements are dependent 
on Mach Number and Reynolds Number. Wire trips 
of different wire diameters were employed in most 
cases in order to obtain transition just upstream of 
separation as flow conditions were varied. In one case, 
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M, = 0.7 and R; = 3 X 108, a 1/16-in.-wide cellophane- 
tape trip was used. The same forebody and trips were 
used for the separated- and attached-flow models. 
Again, tests were performed under nearly identical wind- 
tunnel test conditions for both separated and attached 
boundary layers. 

In contrast to the results for laminar boundary layers, 
the analysis of reference | indicated a large increase in 
heat transfer for turbulent boundary layers when sep- 
arated at subsonic and low supersonic Mach Numbers. 
The analysis also indicated a large effect of Mach Num- 
ber on hs/h4. The results of reference 1 and the ex- 
perimental data are presented in Fig. 7. The experi- 
ments were performed at a constant Reynolds Number 
of about 3 X 10° In contrast to the excellent agree- 
ment of experiment with theory observed for laminar 
boundary layers, a large disagreement is evident for 
turbulent boundary layers, especially at the lower 
Mach Numbers. 

In an attempt to shed some light on this large dis- 
crepancy, velocity and temperature profiles were 
measured in the separated, turbulent boundary layers 
with both heated and adiabatic walls. A shadowgraph 
taken during a probe survey at J/, = 2.6 is presented 
in Fig. 8 along with the pertinent thermal data. These 
data are typical of the probe studies conducted at other 
Mach Numbers and Reynolds Numbers. The tem- 
perature data are presented as the difference between 
the total temperatures measured by the probe with the 
wall heated and adiabatic. The plot in Fig. 8 shows 
how this temperature difference varied through the 


1959 


separated boundary layer. As discussed previously, 
reference 1 assumed that air in the low-velocity, reverse. 
flow region was drawn into the separated boundary 
layer at a temperature equal to wall temperature. [If 
this assumption of the theory were correct, the probe 
would have measured 7, = 7), at the zero velocity 
stream line. The lower portion of the probe survey 
indicates that the measured temperature difference is 
an order of magnitude less than that assumed by the 
theory. This large discrepancy in the temperature 
potential is one possible explanation of the disagree- 
ment found for the experimental and theoretical heat- 
transfer results. The probe surveys also show that 
the thickness of the separated boundary layer is greater 
than one-half the depth of the separated region, which 
may contribute to the difference since it is assumed 
in reference 1 that the thickness is relatively small. 
In view of these detailed studies, it is believed that 
theory must be modified to take into account the details 
of the reverse-flow region. 


Effect of Reynolds Number 


The effect of Reynolds Number on heat transfer in 
separated boundary layers is presented in Fig. 9. The 
data are plotted as average Stanton Number versus 
Reynolds Number for laminar and turbulent boundary 
layers. The results for laminar separated and attached 
boundary layers exhibit the same dependence on Reyn- 
olds Number, St « R,~!/2, which is a well-known result 
for attached laminar boundary layers. Thus, the ef- 
fect of separation on heat transfer is apparently inde- 
pendent of Reynolds Number, as discussed previously 
in this paper. The data for the attached turbulent 
boundary layer exhibit the usual Reynolds Number 
dependence, St, « R,~'/*. The separated-flow data, 
however, exhibit the trend, Sts « R,~*, indicating that 
the reduction in heat transfer by separation increases 
with increasing Reynolds Number. 

In order to explore further the dependence of Stanton 
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Number on Reynolds Number, three additional tests 
were petformed. The original configuration of Fig. 9 
was tested at M/, = 0.3, and a larger diameter unheated 
forebody and afterbody were fitted to the two original 
heated centerbodies and data were obtained at 7, = 
26. The resulting model configurations and experi- 
mental data are presented in Fig. 10. In each case, 
the Stanton Number exhibits a variation with Reynolds 
Number, Sts « R,~?/° in agreement with the earlier 
result. 

A difference of about 30 per cent in the level of Stan- 
ton Numbers for the lower model, compared to the 
upper and center models, is apparent in Fig. 10 (data 
for 17, = 2.6), which requires some comment. As 
mentioned earlier, the lower model was assembled by 
adding a larger diameter, unheated forebody and after- 
body to the heated centerbody of the original separated- 
flow model (top of Fig. 10). As a result, the lower 
model differed from the upper model in two respects: 
(a) the ratio of the depth of the separated region to the 
diameter of the separated boundary layer, Y/D, was 
considerably larger and would presumably result in 
lower reverse-flow velocities; (b) although not apparent 
in the sketches, portions of the lower model surface 
enclosed by separated flow (annular areas under re- 
attachment and separation) were not heated whereas 
the entire surface between separation and reattachment 
was heated for the upper model. Both of these differ- 
ences would probably result in lower measured Stanton 
Numbers for the lower model and the data, in fact, 
show considerably lower values for the lower model. 
The center model also differed from the upper model 
in two respects: (a) the ratio Y/D was less in this case, 
which would presumably result in higher reverse-flow 
velocities and, therefore, in higher Stanton Numbers; 
and (b) since the same unheated forebody and afterbody 
were used on the center model as were used on the lower 
model, portions of the model surface between reattach- 
ment and separation were unheated and would pre- 
sumably tend to decrease the measured Stanton Num- 
bers. It is possible that the opposing effects of the 
two differences, in this particular case, nearly compen- 
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Fic. 11. Streamwise distribution of heat transfer. Laminar 
flow M, = 4, Rt = 0.4 X 108. 


sate since the data of Fig. 10 indicate approximately 
equal Stanton Numbers for the center and upper models. 


Streamwise Distribution of Heat Transfer 


Knowledge of the distribution of heat-transfer rates 
over a surface is important to proper vehicle design. 
All of the data presented thus far have represented the 
average heat-transfer coefficients for constant wall 
temperature. The heated, axially symmetric center- 
bodies were equipped with eight heated segments of 
equal surface area. The measured distributions of 
heat to these eight segments were found to be repeatable 
and are believed to offer a good qualitative indication 
of the variation of local heat transfer on the heated 
surfaces. 

Fig. 11 presents measured distributions of heat 
transfer for laminar separated boundary layers on the 
three configurations shown in Fig. 10. Data obtained 
at 7, = 4.0 for separated flow are presented as the 
ratio of the local heat transfer to the average over the 
heated surface, gs/gs. For reference, the heat-transfer 
distribution for the attached boundary layer, 94/Qa, 
tested at the same outer-edge flow conditions, is shown 
as the dashed line. For the attached boundary layer, 
the maximum heat transfer was tound to occur at the 
upstream portion of the heated surface and the local 
heat transfer decreased in the downstream direction 
which is the usual variation for attached boundary 
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layers. For each of the three models with separated 


boundary layers, however, the maximum heat transfer 
was found to occur at the portion of the heated surface 
farthest downstream. This result of maximum heating 
in the reattachment region was also observed at all 
Mach Numbers and Reynolds Numbers of the present 
tests. 

The streamwise distribution of heat transfer was 
also measured for turbulent separated boundary layers 
on the same configurations discussed above. The re- 
sults for turbulent boundary layers are presented in 
Fig. 12. Again, maximum heat transfer for the 
attached boundary layer is observed at the upstream 
portion of the heated centerbody. Separation of the 
turbulent boundary layer also resulted in maximum 
heat transfer in the region near reattachment for each 
of the three configurations, as shown in Fig. 12. 

Some variations in the shape of the gs/@s plots are 
evident for the various model geometries and the de- 
tailed distribution of heat transfer may be influenced by 
the geometry of the separated region. The result of 
maximum heat transfer occurring in the reattachment 
region was observed for all separated flows of the 
present tests. 


Transition of Separated Boundary Layers With 
Wall Temperatures Below Adiabatic 


Practical applications of boundary-layer separation 
as a means of reducing heat transfer in laminar flow re- 
quire a knowledge of the transition characteristics of 
separated boundary layers. This portion of the paper 
describes some measured effects of Mach Number and 
wall temperature on the Reynolds Number of transition 
of separated boundary layers. For the purposes of the 
present discussion, the transition Reynolds Number 
is defined as the maximum Reynolds Number, based 
on the length of separated boundary layer, ¢’, for which 
transition is just downstream of the reattachment 
region. 

A strong influence of Mach Number on transition 
Reynolds Number was reported in reference 7, for the 
special case of adiabatic wall conditions. These data 
are presented in Fig. 13 as the shaded area. This 
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area encompasses data obtained on a variety of two- 
dimensional models of which the sketches in Fig. 13 
are typical. An increase in Mach Number is observed 
to increase the transition Reynolds Number from the 
order of 50,000 at 7 = 0.3 to about 400,000 at M7 = 3. 
This favorable effect of Mach Number may make sep- 
arated laminar boundary layers of practical interest, 
especially at hypersonic Mach Numbers. 

A knowledge of the effect of wall cooling, 7),/7\,, < 
1.0, on transition is required since aerodynamic heating 
at hypersonic Mach Numbers is usually associated with 
wall temperatures much lower than the adiabatic wall 
temperature. In order to obtain this information, an 
experimental investigation was undertaken recently at 
Ames Research Center using ogive-cylinder models 
cooled internally with liquid nitrogen to a constant 
wall temperature. Some of the resulting data ob- 
tained by Stephen J. Keating, Jr., and the present writer 
are included in Fig. 13. Two separation lengths were 
required because of wind-tunnel stagnation pressure 
limitations. The data obtained for constant values of 
T/Ta,e exhibit a favorable effect of Mach Number 
similar to that reported in reference 7. The effect of 
wall cooling at a constant Mach Number, however, is 
seen in Fig. 13 to decrease the transition Reynolds 
Number for separated boundary layers. This de- 
stabilizing effect of wall cooling is opposite to the sta- 
bilizing effect on attached boundary layers observed 
for moderate wall cooling by other investigators— 
e.g., references 8 and 9. The variable 7), data pre- 
sented in Fig. 13 represent only a fraction of the tem- 
peratures and models investigated. The stabilizing 
effect of Mach Number and the destabilizing effect of 
T. < Tey were observed in all cases. 


Results of Other Investigations 


Before the results of the various experiments (refer- 
ences 2-6) are compared, the types of flows studied 
will be briefly described. Rabinowicz? measured local 
heat-transfer rates on the base of a spike-hemisphere- 
cylinder for turbulent flow in a shock tube. Bloom 
and Pallone* measured local heat-transfer rates on two 
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afterbodies downstream of a blunt nose in a shock 
tunnel, and, for that investigation, transition appar- 
ently occurred in the separated boundary layer. 
Seban* and Emery® measured local heat-transfer rates 
in turbulent separated and reattached subsonic flows 
on two-dimensional plates downstream of a rod and two 
rearward-facing steps. Powers, Stetson, and Adams® 
measured local heat-transfer rates on the base and 
sting downstream of a hemisphere-cylinder using a 
shock tube. The separated boundary layers studied 
in reference 6 were either turbulent or laminar with 
transition apparently occurring in the reattachment 
region. 


The experimental results of the above described 
studies are generally in agreement with the present 
paper, considering the variety of flows investigated. 
A comparison of the results of references 3-6 with the 
analysis of reference | is not discussed since the total 
or average heat transfer between separation and reat- 
tachment could not be calculated from their results. 
In each of these four investigations, approximately 
twofold variations in local heat transfer were measured 
in the regions of separated flow and are in qualitative 
agreement with the present data. In references 2 and 
3, the separated boundary layers did not reattach to 
model surfaces. The two-dimensional separated turbu- 
lent boundary layers investigated in references 4 and 5 
reattached on flat plates and the maximum heat- 
transfer rates were observed in the reattachment re- 
gions, in agreement with the results of the present paper. 
Reference 6 measured no distinguishable variations in 
local heat transfer within the separated region for either 
turbulent or transitional boundary layers. The meas- 
ured average heat transfer for transitional flow was 
compared, in reference 6, with theoretical calculations 
of average heat transfer for attached flow reduced ia 
accordance with reference 1. The comparison shows 
reasonable agreement between the data and theory 
considering the scatter of the data and the probable 
presence of transition in the reattachment region. 
The average heat-transfer measurements of reference 
6 for separated turbulent flows were found to be about 
one-third of the theoretical calculations of that reference 
for equivalent attached turbulent boundary layers. 
This 67 per cent reduction in heat transfer due to 
separation is somewhat larger than the 50 per cent re- 
duction reported in the present paper. This discrep- 
ancy may be the result of the considerably higher stag- 
nation enthalpy levels of reference 6 and differences in 
the attached-flow heat-transfer values used for com- 
parison (measured values in the present paper and theo- 
retical values in reference 6). 


The experiments of reference 6 gave no indication of 
high local heat-transfer rates in the reattachment region 
whereas relatively high local rates were observed for 
the present tests, and in references 4 and 5. One pos- 


sible explanation may be the different types of axially 
symmetric separated flows studied in reference 6 and 


in the present paper. 


The conservation of boundary- 
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layer mass would require a relatively thicker boundary 
layer at reattachment for the converging axially sym- 
metric flow of reference 6 compared with the parallel 
flow of the present experiments, and the two-dimen- 
sional tests of references 4 and 5. 


Concluding Remarks 


Average heat transfer for separated laminar boundary 
layers was found to be about 56 per cent of that for 
eguivalent attached laminar boundary layers, inde- 
pendent of Mach Number and Reynolds Number, in 
agreement with the analysis of reference 1. Average 
heat transfer for separated turbulent boundary layers 
was found to be about 60 per cent of that for equivalent 
attached boundary layers at a Reynolds Number of 
3 X 10° and over a Mach Number range of 0.3 to 4.0— 
in considerable disagreement with the calculations of 
reference 1. The reduction in heat transfer due to 
separation of turbulent boundary layers was found to 
increase with increasing Reynolds Number within the 
limits of the present investigation and, if this trend 
continues to higher Reynolds Numbers, considerable 
reductions in average heat transfer may be possible 
through the mechanism of turbulent boundary-layer 
separation. 


Data are also presented which demonstrate a favor- 
able effect of increasing Mach Number and an adverse 
effect of decreasing wall temperature on separated 
boundary-layer transition. From these observations, 
the practicality of using separation as a means of re- 
ducing heat transfer for laminar boundary layers will 
probably depend strongly on the Mach Number of the 
flow and the ratio of wall temperature to adiabatic wall 
temperature. 


In many applications, it may be desirable to have 
some control of the distribution of local heat-transfer 
coefficients. In all cases of the present investigation, 
maximum heat transfer was measured in the down- 
stream, or reattachment, region of separated flows, 
whereas for attached flows, maximum heat transfer 
was measured at the farthest upstream portion of the 
heated surface as expected. In view of these results, 
single or multiple boundary-layer separations over a 
surface may permit some control of the distribution of 
local heat transfer and result in a reduction of the 
average heat transfer as well. 


References 


1 Chapman, Dean R., A Theoretical Analysis of Heat Transfer 
in Regions of Separated Flow, NACA TN 3792, October, 1956. 

2 Rabinowicz, Josef, Measurement of Turbulent Heat Transfer 
Rates on an Aft Portion and Blunt Base of a Hemisphere-Cylinder 
in the Shock Tube, GALCIT Hypersonic Research Project Memo 
No. 41, 1957. 

3 Bloom, Martin H., and Pallone, Adrian, Shroud Tests of Pres- 
sure and Heat Transfer Over Short Afterbodies With Separated 
Wakes, Dept. of Aero. Eng. and Appl. Mech., Polytechnic Insti- 
tute of Brooklyn (also available as WADC Tech. Note &8-185, 


1958). 


ai 
fer- 
died 
i 
ocal 
ere- 
q 
two 
| 
ers. 


738 JOURNAL OF THE AERO/SPACE 


4Seban, R. A., Heat Transfer and Flow with Separated and 
Reattached Boundary Layers as Produced by Surface Irregularities, 
WADC Tech. Rep. 56-217, 1956. 


5 Emery, A., Heat Transfer Downstream of a Step in a Plate, 
M.S. thesis, University of California, Berkeley, Calif., 1958. 


6 Powers, William E., Stetson, Kenneth F., and Adams, Mac C., 
A Shock Tube Investigation of Heat Transfer in the Wake of a 
Hemisphere-Cylinder, with Application to Hypersonic Flight, Res. 
Rep. 30, Avco Research Laboratory, 1958. 


SCIENCES—NOVEMBER, 1959 


7 Chapman, Dean R., Kuehn, Donald M., and Larson, Howard 
K., Investigation of Separated Flows in Supersonic and Subsonic 
Streams With Emphasis on the Effect of Transition, NACA TN 
3869, March, 1957. 

8 Jack, John R., Wisniewski, Richard J., and Diaconis, N. §, 
Effects of Extreme Surface Cooling on Boundary-Layer Transition, 
NACA TN 4094, October, 1957. 

® Diaconis, N. S., Wisniewski, Richard J., and Jack, John R, 
Heat Transfer and Boundary-Layer Transition on Two Blunt 
Bodies at Mach Number 3.12, NACA TN 4099, October, 1957, 


An Experimental Study of Jet-Flap Compressor Blades 


(Continued from page 702) 


2 Falk, Theodore J., Rotating Stall in Single-Stage Axial Flow, 
Compressors, M.S. thesis, Graduate School of Aeronautical 
Engineering, Cornell University, Ithaca, N.Y., September, 1956. 

3 Sinnette, John T., Jr., and Costello, George R., Possible 
Application of Blade Boundary-Layer Control to Improvement 
of Design and Off-Design’ Performance of Axial-Flow Turboma- 
chines, NACA TN 2371, May, 1951. 

4 Katzoff, S., and Hannah, Margery E., Further Comparisons 
of Theoretical and Experimental Lift and Pressure Distributions 
on Airfoils in Cascade at Low-Subsonic Speed, NACA TN 2391, 
August, 1951. 

5 Sheets, H. E., The Slotted-Blade Axial-Flow Blower, ASME 
Diamond Jubilee Annual Meeting, Paper No. 55-A-156, Novem- 
ber, 1955. j 

6 Lockwood, Vernard E., Turner, Thomas R., and Riebe, 
John M., Wind Tunnel Investigation of Jet-Augmented Flaps on a 
Rectangular Wing to High Momentum Coefficients, NACA TN 
3865, December, 1956. 

7 Clark, Edward L., Jr., An Experimental Study of Jet-Flap 
Compressor Blades, M.S. thesis, Graduate School of Aeronautical 


Engineering, Cornell University, Ithaca, N.Y., February, 1958, 

8 Cassidy, J. F., Design and Calibration of a Cascade Wind 
Tunnel, Graduate School of Aeronautical Engineering, Cornell 
University, Ithaca, N.Y. (to be published). 

® Herrig, L. Joseph, Emery, James C., and Erwin, John R., 
Effect of Section Thickness and Trailing-Edge Radius on the Per- 
formance of NACA 65-Series Compressor Blades in Cascade at 
Low Speeds, NACA RM L51J16, December, 1951. 

1 Herrig, L. Joseph, Emery, James C., and Erwin, John R,, 
Systematic Two-Dimensional Cascade Tests of NACA 65-Series 
Compressor Blades at Low Speeds, NACA RM L51G31, September 
1951. 

1 Erwin, John R., and Emery, James C., Effect of Tunnel Con- 
figuration and Testing Technique on Cascade Performance, NACA 
Rep. 1016, 1951 (formerly, NACA TN 2928). 

12 Dean, R. C. Jr., Aerodynamic Measurements; Gas Turbine 
Laboratory, Massachusetts Institute of Technology, 1953. 

18 Scholz, N., Uber die Durchftihrung systematischer Messungen 
an ebenen Schaufelgittern, ZFW, vol. 4, No. 10, October, 1956. 


The Equation of Energy Balance for Fluttering Systems. . . 


(Continued from page 716) 


16 Lanczos, C., The Variationa! Principles of Mechanics; Univer- 
sity of Toronto Press, Toronto, 1949. 

17 Neumark, S., Concept of Complex Stiffness Applied to Prob- 
lems of Oscillations with Viscous and Hysteretic Damping, Brit. 
RAE Rept. No. Aero 2892. 

18 Ashley, H., and Zartarian, G., Piston Theory—A New Aero- 
dynamic Tool for the Aeroelastician, Journal of the Aeronautical 
Sciences, Vol. 23, No. 12, pp. 1109-1118, December, 1956. 

19 Dugundji, J., Effect of Quasi-Steady Air Forces on Incom- 
pressible Bending-Torsion Flutter, Journal of the Aeronautical 
Sciences, Vol. 25, No. 2, pp. 119-121, February, 1958. 

2” Acum, W. E. A., Aerodynamic Forces on Rectangular Wings 
Oscillating in a Supersonic Airstream, Brit. ARC 13,286, August, 
1950. 


21 Walsh, J., Zartarian, G., and Voss, H. M., Generalized Aero- 
dynamic Forces on the Delta Wing with Supersonic Leading Edges, 
Journal of the Aeronautical Sciences, Vol. 21, No. 11, pp. 739-748, 
November, 1954. 

22 Bisplinghoff, R. L., Ashley, H., and Halfman, R. L., Aero- 
elasticity; Addison-Wesley Publ. Co., Inc., Cambridge, Mass., 
1955. 

23 Pines, S., An Elementary Explanation of the Flutter Mecha- 
nism, Proc. Nat. Special. Meeting, IAS, Forth Worth, Tex, 
November 6-7, 1958. 

24 Pulos, John G., An Energy Solution to the Three-Dimensional 


Flutter Problem, thesis for degree of Master of Science in Applied 
Mechanics, Polytechnic Institute of Brooklyn, June, 1949. 


earlier ). 


Microcards for Volumes 1-5 Available 


Microcards for Volumes 3-5 of the JOURNAL are now available (Vols. 1 and 2 were republished in this form 
All five volumes are priced at $5.59 each, and are obtainable through: 


J. S. Canner & Company, Inc. 
Boston, Roxbury 20, Mass. 


The 
trary 
directi 
proxin 
circun 
influet 
of dist 
influer 
illustré 
resent 
ing inl 
signific 

analys 
induce 
a trav 
cludin; 
ever, 1 
purely 
hood, 
cipient 
tremel 
eters a 


Rece 

Th 
ratory 
summe 
R57AG 
his gra 
cussion 
Klappr 
throug] 

* Pre 


53 
Ait 
4 
fa 
by 
4 
L, 
m 
M 
n 
N 
N* 
> 
I 
p 
H 
H 
2) 


Howard 
subsonic 


CA TN 


» N. 
ition, 


ohn R,, 
Blunt 
957, 


, 1958, 
Wind 
~ornell 


hn R,, 
te Per. 
cade at 
hn R., 
-Series 


ember 


Con- 
NACA 


urbine 


ungen 
56. 


Aero- 
Zdges, 
+-748, 


Aero- 
fass., 


‘echa- 


ional 


plied 


Summary 


The flow behind a rotor (or stator) blade row due to an arbi- 
trary inlet distortion along both circumferential and spanwise 
directions is treated on the basis of classical actuator dise ap- 
proximation. It is shown that, in general, the spanwise and the 
circumferential distortion waves interact, so that the combined 
influence is not the sum of the influences of the two simple types 
Numerical examples are then given to show the 
Also 


of distortions. 
influence of the various conventional design parameters. 
illustrated is the role played by the cascade characteristics, rep- 
resented by the variations of loss and exit flow angle due to chang- 
ing inlet flow angles. It is shown that both quantities exert a 
significant influence on the flow under inlet distortion. 

Using the same approach to the problem of incipient stall, the 
analysis shows that theoretically there can be two types of self- 
induced distortion—a purely spanwise, axisymmetric type, and 
a traveling, circumferential plus spanwise type, the latter in- 
cluding a purely circumferential tvpe as a special case. How- 
ever, the analysis also indicates the interesting result that the 
purely circumferential type of traveling stall will, in all likeli- 
hood, precede the other types. One aspect of delaying the in- 
cipient traveling (rotating) stall can then be reduced to two ex- 
tremely simple rules with respect to the desirable design param- 
eters and cascade geometry. 


Symbols 

a = tana 

b = tanB=k+a 

k = a ratio of velocities. (With inlet distortion kU is 
the peripheral speed of rotor. With self-induced 
distortion kU is the speed of distortion propaga- 
tion. ) 

: a = circumferential length of cascade, hence also pri- 
mary wave length of circumferential distortion 

L. = cascade blade span 

m = number of half distortion waves in spanwise direc- 
tion 

M = defined by Aa = 2M(Ap;) 

n = number of full distortion waves in circumferential 
direction 

N = defined by A tan 8. = NA tan £6, 

N* = defined by Ap, = N*A8, 

FP = pressure associated with the main flow 

p = pressure perturbation 
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P, = total pressure associated with the main relative 
flow 

pt = total-pressure distortion or perturbation 

R = reaction ratio 

R = recovery ratio = 1 — (pp/pu) 

ee = main absolute velocities (see Fig. 1) 

AU, AV = absolute flow distortion 

u,v,Ww = induced velocities 

x, ¥,2 = coordinates (see Fig. 1) 

a), a2 = absolute main flow angles, upstream or downstream 
of cascade 

By, Be = relative main flow angles, upstream or downstream 
of cascade 

AB, AB2 = change of relative flow angles immediately up- 
stream or downstream of cascade 

6 = time lag in the response of Aé and Ap» with As, 

€ = phase change of exit distortion compared to inlet 
distortion 

An = wave number for circumferential distortion = 
2an/Ly 

Hive = wave number for spanwise distortion = mx/L, 

® = total-pressure loss coefficient of cascade 

¢ = velocity potential of induced velocities 

= flow coefficient 

v = work coefficient 

Subscripts 


1 for upstream of cascade 
2 for downstream of cascade 
3 for downstream of a stator in a complete rotor-stator stage 


Introduction 


a the aerodynamic design of compressors 
today is based mainly on an axisymmetric and 
steady flow model, there are certain nonaxisymmetric 
and nonsteady flow phenomena which can no longer be 
ignored. Among these are inlet flow distortion and ro- 
tating stall. 

One type of inlet distortion can still be considered 
axisymmetric and steady; this is the purely spanwise 
distortion analyzed by Smith.' Another type, that of 
a purely circumferential distortion,** has been investi- 
gated by Ehrich? and by Rannie and Marble.’ It 
appears that a unified treatment of distortions consist- 
ing of both spanwise and circumferential wave patterns 


** The author’s attention was called by the reviewer of this 
paper to a recent publication on purely circumferential inlet 
distortion by Katz!! who presented some interesting analytical 
and experimental results for the case of a rectangular depression 
of the upstream total pressure. 
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Fic. 1. Isolated plane cascade of finite span with spanwise and 


circumferential inlet distortions. 


is possible; such a treatment-would then include purely 
spanwise and purely circumferential distortions as 
special cases. It is shown in this analysis that the 
spanwise and the circumferential distortion waves 
interact so that the combined influence is not the sum 
of the influences of the two simple distortions. Nu- 
merical examples are shown for various values of the 
conventional design parameters—i.e., flow coefficient, 
work coefficient, and reaction ratio. Also illustrated 
is the influence of the cascade characteristics, repre- 
sented by the variations of loss and exit flow angle with 
a changing inlet flow angle. 

Compared to inlet distortion, the problem of rotating 
stall has received much greater attention since its first 
description by Iura and Rannie.* Various theories 
have been proposed, the majority of which use a small 
perturbation approach,* such as the ones by Emmons, 
et al.,5 Sears,® Marble,’ Stenning, et al.,*> Yeh,’ and 
Rannie and Marble.* Since a small perturbation ap- 
proximation is used in the theories mentioned above, it 
is clear that they apply strictly for an incipient stall 
which is identified in reference 3 as a self-induced dis- 
tortion. 

It is shown in this analysis that theoretically there 
can be two types of self-induced distortion—a purely 
spanwise and, therefore, axisymmetric type and a 
traveling, circumferential plus spanwise type, the latter 
including a traveling, purely circumferential type as a 
special case. It is further shown that the purely cir- 
cumferential type would in all probability occur first. 
The influence of the conventional design parameters 
on self-induced distortion is then discussed in some 
detail. 


* Analyses without the assumption of small perturbation have 
been made by Fabri and Siestrunck'? and Kriebel.'* See a liter- 
ature survey by Emmons, Kronauer, and Rockett.' 
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Analysis 


We shall consider the case of a plane cascade of in- 
finite blades and of finite span between two walls (Fig. 
1). The cascade itself is approximated by an actuator 
disc located at x = 0. The coordinate axes are chosen 
as follows (see Fig. 1): x represents the axial direction, 
y the circumferential direction, and z the spanwise direc- 
tion. The region upstream of the cascade is denoted 
by the subscript 1, and that downstream of the cascade 
by subscript 2. (Equations without subscripts apply 
to both regions.) 

The undistorted or main absolute velocities are de- 
noted by U along x and V along y. The absolute main 
flow angle a is, therefore, a = tan~! V/U. The rotor 
blade velocity is represented by kU, where is the re- 
ciprocal of the conventional flow-coefficient ¢. (For 
a stator, k = ().) The relative main flow angle @ is, 
therefore, 8 = tan-'[(V + kU)/U]. It is clear that 
tan @ = k+ tana. The main flow is assumed to be 
uniform along the span (as well as in the circumferential 
direction), hence the blades are untwisted. However, 
the turning of main flow across the cascade need not be 
small. 

The inlet distortion is represented by the velocity 
components AU; along x and AV; along y far upstream 
atx = —o. These velocity components are, in gen- 
eral, functions of both y and z; in a turbomachine this 
means a combined circumferential and spanwise inlet 
distortion. Due to this inlet distortion, the velocity 
far downstream at x = + will be distorted by an 
amount AU, and AV. (At present we shall leave the 
question open as to whether AV’, = AU) tan a, Al: = 
AU, tan a.) Thus; 


Velocity far upstream: 


U+ AU(—@, 9,2), Vit y, 2), 


Velocity far downstream: 


U + y,2), Vot+ AV2(+, y, 

In the following it will be assumed that the velocity 
distortion is small compared to the main velocities. It 
is then possible to consider that, in the entire upstream 
or downstream region, the velocity perturbation con- 
sists of (1) the distortion velocities AU’ and AV carried 
along by the main velocities such that 


AU = AU(x tan a — y,z), AV = AV(xtana — y, 2) 
and (2) an additional induced velocity field represented 
by 

u= u(x, y, 2), v= v(x, 2), w= w(x, y, 2) 
The purpose of this analysis is to find this induced ve- 
locity field which will then yield the exit distortion 
AU, and AV», all under a given inlet distortion AU, and 


AV, which, as will be seen later, may be arbitrarily 
prescribed. 
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Basic Equations 
Substituting the velocity components 


U + AU(x tan a — y, 2) + u(x, y, 2) 

V + AV(x tan a — y, 2) + v(x, ¥, 2) 

w(x, 2) 
into the incompressible, nonviscous equations of motion 
and retaining only the first-order terms, we obtain 


U(Ou/Ox) + V(Ou/Oy) = —(1/p)(Op/Ox) (1) 
U(dv/Ox) + V(dv/dy) = —(1/p)(Op/dy) 
U(Ow/dx) + V(Ow/dy) = —(1/p)(Op/dz) (3) 


In the above equations, p = p(x, y, Z) is the pressure 
perturbation. 

The continuity equation for incompressible flows 
gives 


(Ou/Ox) + (Ov/Oy) + (Ow/dz) = 0 (4) 


After canceling out p in any pair of the first three 
equations, we obtain three equations of the form 


U[o( )/ox] + )/oy] = 0 


where ( ) represents any one of the three vorticity 
components 


[(Qu/Oy) — (Ov/Ox)], [(Ov/Oz) — (Ow/Oy)], 


[(Ow/Ox) — (Ou/dz) 


Physically it merely means that the vorticity is car- 
ried by the main velocity only and not by the induced 
velocities, a well-known result of linearization. Since 
u,v, ware zero atx = +, their vorticity must be zero 
everywhere. Hence, we can define a potential ¢ such 
that 


u = 0¢/dx, v = O0¢/dy, w= I¢/dz (5) 


Substituting into the continuity equation, the result 
is a Laplace equation for g. Keeping the appropriate 
boundary conditions in mind, the solution for ¢ up- 
stream of the cascade and for ¢g. downstream of the 
cascade are as follows: 


m n 2%. (6) 
= DY Vem? + + Any COS 
m n 
where 


blade span, m = 0, 1, 2,3,... 


Hm = L; 


primary wave length along y, 
0, 


A, = 2nx/L,, Ly, 


The coefficients 1,,, and B,,, are in general complex 
numbers. From Eq. (5), 


Amn V sim? + rn? 4 


m n 


uy 


eV um? + drn?x + COS (7) 


a= > A mntdneV + An?x + COS (8) 
m 


n 


w= — > A mnblme dnt x + sin (9) 


m n 


with similar expressions for “2, v2, We. It is seen that 
w = () at the walls (where z = 0, L,) as required. 


Matching Conditions at Actuator Disc 


For the remaining part of this analysis, we shall use 
the subscript 1 for x = —0 and subscript 2 for x = 
+0, these being the stations immediately upstream and 
immediately downstream of the actuator disc. Con- 
sidering a single upstream two-dimensional wave in y 
and z represented by 


Any 
1 = AmnV bm? + An? COS (10a) 


with corresponding expressions for 7, w;, and a single 
downstream two-dimensional wave represented by 


= BV + COS (10b) 


with corresponding expressions for v2, we, we obtain, 
from Eqs. (7) and (8), the following important relations: 


= 1Vmnlti, V2 = (11) 


where 


= bm? 4 pore 
(2n/Ly)/V (2n/Ly)? + (m/L,)? 


The first matching condition at the disc is continuity. 
This requires 


AU, + = AU U2 (12) 


The second matching condition is that the spanwise 
velocity should undergo no change across the actuator 
disc—i.e., 


W, = We (13) 


This condition is probably unimportant when yu, < A,, 
since w would then be much smaller than the other in- 
duced velocities u and v. (Notice that w is identically 
zero at the walls.) On this aspect more will be said 
later. 

If Eq. (13) is valid, then a comparison of Eq. (9) for 
w, and a corresponding equation for w. shows that 


Aun = tims. Then, 
= = Ve (13’) 
For m = 0), Eq. (13) is identically satisfied since w, = 


We = 0. 


Cascade Characteristics 


In steady-state cascade testing, a cascade of fixed 
geometry (stagger, solidity, camber, etc.) is subjected 
to various inlet flow angles 6;. Considering 6, as the 
independent variable, the following quantities are the 
dependent variables: (1) exit flow angle 62, (2) loss of 
total pressure (this is usually represented by a loss 
coefficient @ defined as the ratio of total-pressure loss 
to inlet dynamic pressure), and (3) static pressure rise. 
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Any two of the above three quantities are sufficient to 
determine the third. For our purpose we shall choose 
the first two—i.e., and 

Consider first the exit flow angle 62. For small vari- 
ations, we can represent a function by a linear relation. 


Thus, 
AB, = (14) 
where N* is the slope in the 6, versus 6; curve at the 


cascade operating point. It is fairly simple to relate 
AB, and AB, to the velocity perturbations. From 


Ut AV +1 


we get 
A tan = [((AV + v)/U] — [((AU + u)/U] tan B 
Hence, 


AB = cos? B(A tan B) = 
(cos? B/U) [((AV + v) — (AU + x) tan B] 


Eq. (14) then gives 


(AV2 + %) — (AU2 + uz) tan = 
N[(AVi + a1) — + m) tan (15) 


where NV = N* cos? 6,/cos? Notice that Eq. (15) 


is equivalent to 
A tan = N A tan 


Hence, NV is the slope in the tan 62 versus tan 6; curve. 
At design condition JN is usually small, but may not be 
negligible especially for low solidities. 

Making use of Eq. (11), Eq. (14) or (15) can finally be 
written as 
AU; tan Bo + (tan By + 1Vmn) U2 

N[AU, tan By AV; + (tan Bi 1Vmn)U1 | (16) 

Next consider the change of relative total pressure 
across the cascade. (By relative total pressure we 
mean the sum of the static pressure and the dynamic 
pressure due to velocities relative to the blade.) The 
relative total pressure at either side of the cascade con- 
sists of a main relative total pressure P, (due to main 
flow alone) and a total-pressure perturbation p,. Thus, 


Pit pi= P+ p+ (1/2)e((U + AU + u)? + 
(V + AV + »)?] 
~ P+ p + (1/2)e(U? + V2) + 
p[U(AU + u) + V(AV + »)] 
Since by definition 
P. = P + (1/2)e(U? + V2) 
it follows 
bi = p + p[U(AU + u) + V(AV + 2)] 


But from Eq. (1) and 0u/O0y = 0v/Ox, it can be shown 
that 


p+ p(Uu + Vv) = 0 


(This merely means that the induced flow is a potential 
flow, and, hence, will not change the total pressure.) 
Thus, the relative total-pressure perturbation due to 
the velocity perturbation is 


bi = p(UAU + VAV) (17) 


The loss of relative total pressure across the cascade 
can be expressed as 


(Pu + pa) — (Pe + Pe) = (@ + Ad)(p/2) X 

+ AU, + m)? + (Vi + AV; + 
where @ is the loss coefficient for the main flow and 
A@ is the change of due to the inlet-angle perturba- 


tion as a result of velocity perturbations. Now, for the 
main flow alone, 


Pu — Pe = &(p/2)(U? + Vi") 
Hence, 


p 


Pu — Pe = Ao — + tan? Bi) X 


[1 + +m) + tan + 1) 
AGU,(1 + tan? 


bo! 


The second term in the above brackets is of the order of 
@/(A@/246,). It turns out in later calculations that, 
when the loss begins to exert a significant influence, 
A@®/2A8; is of the order 1, whereas @ is still one order of 
magnitude smaller than 1. Thus it appears permissible, 
for the sake of reducing the amount of mathematical 
calculation, to drop the second term in the brackets 
and let 


Pa — Pe = + tan? 


We now substitute p, and p,. in the above equation in 
terms of perturbation velocities by Eq. (17). Further- 
more, for small changes of 6; (i.e., small A6,), A& can be 
expressed as a linear function of A;; i.e., 


Ae = 2MA,, (18) 


where .V/ is one-half of the slope of the @ versus §; curve. 
The resulting equation after these substitutions is 


(AU, + AV» tan — (AU; + AV; tan — 
1Vmnk (Us + = M[AU, tan By AD, + 
(tan By 1Vmn) tr | (19) 


Eq. (16) and (19) represent the cascade characteristics. 


Solutions 


In the foregoing we have derived the three Eqs. (12), 
(16), and (19). In addition, there is also Eq. (13’) for 
m ~ (. Given the inlet distortion AU and AV;, the 
main flow angles ai, a2, 8:, 82, and the cascade char- 
acteristics V and M, we should be able to solve for the 
exit distortions AU; and Al. The best way to do 


this is to consider separately the following four cases. 
To simplify future writing we shall introduce the 
symbols 
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a, = tan a, = tana, = = tan Bo 


Case (1)—-Purely Spanwise Distortion 
In this case 2 = A, = Ym = VU. From Eqs. (12) and 
(13’), 
AY, AU. = = —2u, 
Using the above relation to eliminate m, and mw in Eqs. 
(16) and (19), we obtain two equations for the two un- 
knowns AU.and AV»: 
(bs = Nb) = 
(Nb, — bo)AU, — 2NAV, 
(2 Mb,)AU2 Ve = 
(2 + Mb,)AU, + 2(b; — M)AV, 


(20) 


This pair of equations can easily be solved. The re- 
sults can be considerably simplified if we separate the 
two types of inlet distortions: (a) the distortion ve- 
locity parallel to the main flow, and (b) the distortion 
velocity perpendicular to the main flow. Any other 
types of inlet distortion can then be obtained by a com- 
bination of these two basic types. 

(la) The Inlet Distortion Velocity Parallel to the Main 
Flow—For this case we have AV; = AUja;. The solu- 
tion to Eq. (20) then becomes 


RECOVERY RATIO @ 
50% REACTION M=N=O 


EQ. (24) 


FLOW. COEF. 


! 


Recovery of purely spanwise distortion, 50 per cent 


Fic. 2. 
reaction, I = N = 0, (inlet distortion velocity parallel to 
main flow). 


AUs _ [M(2k — by) + 2(1 + asby)| + — by) 


AU, 2 — Mb, + bx(b, — (21) 
_ (bz — Nbi)[M(2k — bi) + 2(1 + aibi)| — (2 — Mbi)(2RN — Nbi — i 


If in addition V = M = 0, Eqs. (21) reduce to 


AU2/AU, = [2(1 + aybi) — bo?]/(2 + 


AU, 2(2 — Mb, + be(b, — 


AV2/AU, = bo(2 + ayhi)/(2 + by?) (22) 


Thus, although the inlet distortion velocity is parallel to the inlet main absolute velocity, the direction of the exit 
distortion velocity is, in general, neither along the main exit absolute velocity nor along the main exit relative ve- 
locity. On account of this, the exit flow (including the distortion) will have slightly different directions at different 


spanwise locations. 


For the simple case represented by Eqs. (22), it is of interest to examine one complete stage in a multistage ma- 


chine. 


The exit distortion of the rotor given by Eqs. (22) becomes the inlet distortion to the following stator if the 


interaction between the two blade rows is neglected. The exit distortion of the stator can be obtained from Eqs. 


(22) by letting k = 0. 
of repeated stage design (i.e., @3 = a) are 


Using the subscript 3 for stator exit, then the exit distortions Al’;, AV’; for a complete stage 


AU, (2 + a,*)(2 + dy?) 23) 
AVs _ ai[4(1 + aibs) — + + aibs)] 
AU, (2 + a1*)(2 + by?) 


Another quantity of interest is the total-pressure per- 
turbation at rotor exit as a ratio of the same quantity 
at inlet. For the simplified case VN = M = 0 repre- 
sented by Eqs. (22), this ratio is 


Pe + debe) — 


Pa (1 + a,*)(2 + 


A more general expression for this ratio, including the 
influence of the parameters NV and M, can be derived in 
The significance of this ratio is as 


(24) 


a similar manner. 


follows: The closer this ratio is to (, the more is the 
tendency of the stage to smooth out any inlet distortion. 
Smith! has defined a “‘recovery ratio” R=1 — (pr/pu). 
If R = 0, a total-pressure distortion upstream of the 
rotor will be fully carried over to the downstream. If 
G = 1, the upstream distortion will be completely 
smoothed out, and the downstream flow will experience 
no distortion in total pressure. If ® > 1, an upstream 
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total-pressure defect (or excess) will become a total- 
pressure excess (or defect) at downstream. Con- 
versely, if R < VU, an upstream velocity distortion will 
become worse downstream. Thus it seems that a re- 
covery ratio in the neighborhood of 1 would be highly 
desirable in an axial flow turbomachine. 

Since Smith’s results [Eqs. (2), (3), and Figs. 5-7 of 
reference 1] are obtained on the basis of somewhat 
different assumptions, it is interesting to compare his 
results with those of the present analysis. For this 
purpose we shall adopt the following conventional 
parameters: 


Flow coefficient % = 1/k 


Work coefficient = Work/[(kU)?/2] = 
(2/k)(a, — a2) = (2/k)(b, — be) 


Reaction ratio R = (bj + be)/2k = 
1+ + a2) /2k | 


These parameters take the place of the angles a, a, 
62 or their tangents do, db), bo. For 50 per cent 
reaction (Fig. 2), 25 per cent reaction (Fig. 3), and 75 
per cent reaction (Fig. 4), the recovery ratio ® is 
plotted against the flow coefficient ¢ for two values of 
v¥: w = 0 and y = 1. In each case, the solid lines 
are calculated from Eq. (24) and the dotted lines are 
based on Eqs. (7), (8), and (9) of reference 1. It can be 
seen that the two different analyses yield qualitatively 
similar results; although, in the range where the recov- 
ery ratio reaches its peak, its value is higher and occurs 
at somewhat lower flow coefficient in this analysis than 
in reference 1. (It is interesting to note that for 25 per 
cent reaction and y = 1, the two analyses give identical 
results.) As pointed out by Smith, his results show 
that (1) the flow coefficient should be around ().5 or 
somewhat less, (2) light loading is generally desirable 
especially at higher flow coefficients, and (3) the re- 
covery ratio is not a strong function of reaction ratio. 
These results are again approximately valid on the 
basis of this analysis.* One could, however, modify 
slightly the first two conclusions as follows: (1) the 
flow coefficient should be around 0.5 or somewhat more, 
and (2) light loading allows somewhat higher flow 
coefficients for full recovery. 

The above conclusions are for N = M = 0—i.e., 
the cascade exit angle and the loss are assumed to be 
independent of the inlet angle. To find out the in- 
fluence of the variation of loss, two cases M = ().1 and 
M = 1.0 are used to calculate the recovery ratio for 50 


* The essential difference between the two analyses is as fol- 
lows. The analysis of reference 1 assumes that the flow pertur- 
bation at blade trailing edge is the same as that at far downstream, 
whereas this analysis assumes that the flow at blade trailing edge 
is the same as that immediately behind the actuator disc. The 
truth of the matter must lie somewhere between these two as- 
sumptions. Presumably the assumption of reference 1 works 
better for blades with low aspect ratio, whereas the assumption 
in this analysis works better for blades with high aspect ratio. 
It is significant, therefore, that both analyses reach roughly the 
same conclusion with respect to the desirable values for flow 
coefficients. 


COMBINED DISTORTION 
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M=0O 
a A: PURELY CIRCUMFEREN- 
TIAL, M=O, N#0 
B: m=!, N=! 
N=2 
FLOW COEF. 
Fic. 7. Ratio of total-pressure perturbation for combined 


circumferential and spanwise distortion, 50 per cent reaction, 
¥ = 1, N =0,L,/L; = 3n. 


per cent reaction, y = 1, using the exit velocity dis- 
tortions as given by the more general Eqs. (21). 
These two values of ./ correspond to an increase of loss 
coefficient @ by 0.0035 and 0.035, respectively, per 
degree increase of inlet flow angle. (The larger figure 
represents, therefore, a considerable stall.) The re- 
sults are compared to the previous case of \J = 0 in 
Fig. 5. It appears that the influence of a mild .\/ is to 
decrease the recovery ratio somewhat, but does not 
change the main characteristics or the conclusions of 
the previous paragraph. On the other hand, with a 
severe \/ the recovery ratio is negative almost all the 
way. It would appear, therefore, that no practical 
recovery is attainable under severe stall. 

The influence of N is investigated next, again for the 
case of 50 per cent reaction and Yy = 1. Here a value 
of N = 0.2 is compared to N = 0. (N = 0.2 corre- 
sponds to a compressor cascade of solidity around 0.6 
at reasonable stagger and incidence.) It is seen in 
Fig. 5 that N tends to reduce recovery ratio for all 
practical values of ¢. The combined influence of .V 
and .M at design condition seem to suggest that for best 
recovery of spanwise distortions, the flow coefficient 
should be in the neighborhood of 0.4 to 0.5. 

Although the results represented by Eqs. (20) to 
(24) are derived for a single wave given by Eg. (10) 
with m = 0,m = 1,2,3..., the fact that these results 
are independent of m means that all wave lengths of 
the radial distortion are attenuated or amplified by the 
same factor. This, plus the fact that there is no phase 


shift in exit distortion (compared to inlet distortion), 
enables us to conclude that the relative profile of any 
radial inlet distortion is preserved at the exit, and only 
its magnitude is modified by a constant factor. 
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ever, as pointed out previously, the exit flow angle will change by slightly different amounts along the span. As 
will be seen later, these conclusions are quite different from the case of a circumferential distortion. 
(1b) The Inlet Distortion Velocity Perpendicular to the Main Flow*—For this case we have AV; = —AU/a;. The 


solution to Eq. (20) then becomes 


AU, b.[N(2 + = | + [M(2 + ayb;) 2k | 
AU, a;(2 Mob, bo(be Nb) | (25 
— Nb:)[M(2 + abs) — 2k] — (2 — Mby)[N(2 + — 
AU, 2a, [2 Mb, + bo(be Nb,) | 
If, in addition, V = M = 0, 
AU2/AU, = (—aib2? — 2k)/as(2 + bp”), AV2/AU, = bo(a, — k)/a;(2 + (26) 


Since pu = O when the inlet distortion velocity is per- 
pendicular to the main flow, it is impossible to define 
an attenuation factor on the basis of total-pressure 
distortions. One alternative is to define a velocity 
attenuation factor @ as the ratio of the exit velocity 
distortion in the direction of the exit main flow to the 
inlet distortion velocity (which is perpendicular to the 
main flow). For the simple case V = M = 0, 


= — [sin a; cos ay/a;(2 + + + a2)] (27) 


This attenuation factor @ is plotted for 50 per cent reac- 
tion blading, with Yy = Oand y = 1. (See Fig. 6.) It 
is seen that here again the flow coefficient should be 
around 0.5 for the least downstream distortion in total 
pressure, corresponding to a vanishing @. 


Case (2)—Purely Circumferential Distortion 


In this case m = ptm = 0, ¥mn = 1. The three avail- 
able Eqs. (12), (16), and (19) contain four unknowns 
AUs, AVo, Ue. For an inlet distortion AV, = 
a,AU,, which is moreover steady in the absolute frame 
of reference, the flow at exit is also steady. A con- 
sideration on continuity and equations of motion at far 
downstream then shows that 

AV2 = aAU2 (28) 
Making use of this relation, Eqs. (12) and (16) can now 
be used to solve for “4; and w%. The results are 


(2 k— by 1)AU, (a2 1)AU2 


Substituting 1“, and uw, into Eq. (19), we can solve for 
AU 2: 


AU2/AU, = (A + 1B)/(C + iD) (30) 
where 


A = (1 + ab + Mk)(b, — Ndi) + 
Mb,(Nk — be) + (k — RN — M) 
B = (1 + ab + Mk)(1 + N) + 
RN(2k — — M) — hbo + — By) 
(1 + deb.) — Nb) — 
M(1 + ash) + — N) 
1 + + k? + — bi) + 
N(1 + deb, — 


D 


For a stator, k = 0 anda = b. Hence, across a stator, 
AU;/AU2 = (As + 1B,)/(C, + (31) 


where 


A, = (1 + a,*)(a3 — Naz) — M(1 + aeas) 
= (1 + a,”?)(1 + N) + — az) 

(1 + a3")(a; — Naz) — M(1 + 

s = (1 + a;?)(1 + N) + M(az — az) 


SAR 


The product of Eqs. (30) and (31) gives the attenuation 
for a complete stage consisting of a rotor and a stator. 
(For a repeated stage in a multistage machine, a; = 


(b2 — Nb:) + i (1 + N) (9) 

uz = AU, — AU, + um For the specially simple case when both V = M = 0), 

= Eg. (30) becomes 
AU, (1 + + + + — kde) (32) 
AU, (1 + + k + + Rk? + ay”) 1 
This equation can be written in an alternate form: AU2/AU; = |AU2/ AU, |e” (33) 

| | 2 2 52) 1/2 
_ {a + + (1+ (34) 

AU) (1 + az”) [4k? + (1 + 
2\(9 2... obs 

(1 + b2°)(2 + arbi + + 2k? — 2(1 + aidi)(1 + (35) 


* Dr. L. H. Smith, Jr., has pointed out to the author that such distortion velocity may arise due to the upstream secondary flow. 
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Fic. 8. Phase shift across rotor for combined circumferential 
and spanwise distortion, 50 per cent reaction, y = 1, N = 0, 


L,/L,z = 3rn. 


Eqs. (34) and (35) have been obtained by Ehrich.” 
The quantity |AU2/ AU, | represents the attenuation 
on the magnitude of the exit distortion of the rotor, 
whereas e, represents a phase shift in the circumferential 
direction. 

Calculations based on these equations will be included 
in the more general case to be described next. It 
might be remarked that Eqs. (34) and (35) show that 
both the attenuation and the phase shift are independ- 
ent of the inlet distortion wave number ),._ When the 
inlet distortion consists of many different wave num- 
bers, the exit distortion will not in general maintain 
the same shape, since a constant phase shift ¢ for all 
wave numbers means that waves of different wave 
numbers will be displaced relative to one another. 


Case (3)—Combined Spanwise and Circumferential 
Distortion With pn<< Xn 


When the spanwise wave number yu,, is much less 
than the circumferential wave number X,, correspond- 
ing to a distortion velocity gradient which is compar- 
atively mild spanwise and comparatively severe circum- 
ferentially, we would expect the induced velocities 
w<u =v. For such cases it is reasonable to assume 
that Eq. (13) can be ignored. We, therefore, have the 
four equations (12), (16), (19), and (28), from which 
the four unknowns AV2 can be solved. The 
results are as follows. 


Across a rotor: 
AU2,/AU, = (A + iB)/(C + iD) (36) 


where 


A = (1 + + Mk)(b2 — + 
Mbi(Nk — bz) + Ymn2(k — RN — M) 


EB = (1 + aid, + Mk)(1 + N) + 
(2k — M)(Nk — be) + (kbe — kNb; — Mb,)} 
(1 + (be Nb) Mayzb; + 


— RN — M) 
D= (1 + deb, — + 
M(az — bi) + 1 + + 
Across a stator: 
AU;/AU2 = (A, + iB) /(C, + iD) (37) 
where 
A, = (1 + — Naz) — M(vmn? + 
Bb, Vmnt (1 + a*)(1 + N) + M(a; 
C, (1 + a3") (a3 Naz) M (inn? + Q2Q3) 


Ds, = (1 + a3*)(1 + N) + M(az — a2)} 


In particular, if = 0, Eq. (37) shows that 


AU;/AU2, = (1 + as?)/(1 + a3") across a stator (38) 


Comparing the above equations with Eqs. (30) and 
(31) previously derived, it can be seen that Eqs. (30) 
and (31) are a special case of Eqs. (36) and (37) with 
= 1. 

The ratio of total-pressure perturbation is 


= (AU2/AU,)((1 + ay?)/(1 + a*)] = 

le across arotor (39) 
= (AU3/AU2)[(L + + = 

| across a stator (40) 
where AU2/AU, and AU;/AU2 are given by Eqs. (36) 
and (37). Asacheck, one notes that for MW = 0, py = 
Py.and e, = 0 across a stator. 
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Fic. 9. Ratio of total-pressure perturbation for combined 
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and spanwise distortion, ¢ = y = 1, N = 0, L,/L. = 37. 


In order to find the influence of »v,,,, let us consider a 
machine of hub-tip radius ratio of 0.5, which gives 
L,/L, = 3m where L, is the circumference at pitch 
radius and L,isthespan. Then 


Von = Vno = 0, = = ().2076, 


voy = 0.1055, vy = 0.3907, ete. (41) 


Fig. 7 presents | Peo/ bul across the rotor calculated 
from Eqs. (39) and (36) for a 50 per cent reaction stage 
with y = 1, N = 0, and two values of VV; Oand 0.1. 
It appears that all curves reach a minimum value at 
around @ = 0.5. Also, the curves form = n = 1 and 
m = 1, n = 2 are lower than that for m = ( for the 
usual values of ¢, the difference being especially pro- 
nounced near @ = 0.5. Thus, it appears that a com- 
bined spanwise and circumferential distortion is atten- 
uated more than a purely circumferential one. Also, 
the curve for v2 is significantly lower than that for vy, 
meaning that more circumferential ‘‘wiggles’ in inlet 
distortion are helpful in recovery. Fig. 7 also shows 
that a small nonzero JM increases somewhat the exit 
distortion. This is similar to a purely spanwise distor- 
tion. 

Fig. 8 presents the phase shift e, across the rotor as 
defined by Eq. (39). It shows that the phase shift 
may vary from slightly negative to —180° or more 
depending strongly on the value of ¢. For @ > 0.7, 
the case of combined spanwise and circumferential 
distortion does not differ much from that of a purely 
circumferential distortion. However, at @ = 0.5 or 
less, the former may give considerably more phase shift 
than the latter. The introduction of 1/7 appears to 
reduce the absolute amount of phase shift. 

In order to examine the influence of reaction ratio 
R, a calculation is made on @ = y = 1, with varying R 
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from 0 to 1. The pressure perturbation | Puo/ pul is 
shown on Fig. 9 for V = 0 and two values of 7; 0 and 
0.1. The phase shift for these cases are shown on Fig. 
10. It appears that the exit distortion increases with 
reaction ratio; lower reaction thus appears preferable 
on this account. The interaction of spanwise distortion 
with circumferential ones shows up in a rather compli- 
cated manner. The influence of V/ is to increase the 
exit distortion. Considerable variation of phase shift 
is exhibited with respect to reaction ratio R and ppp. 

Since a small value of 1/(= 0.1) does not appear to 
exert influence on | and phase shift, a 
larger value of J = 1 is chosen, for the case of purely 
circumferential distortion, 50 per cent reaction, y = 1, 
N= 0. The results are plotted on Fig. 11 for 
and on Fig. 12 for the phase shift. As far as the dis- 
tortion magnitude is concerned, Fig. 11 shows that for 
@ < 0.62, M reduces the exit distortion and hence gives 
better recovery, whereas the reverse is true for ¢ > ().62. 
Fig. 12 shows that the absolute magnitude of phase 
shift is decreased considerably by a large M. 

The influence of the variation of exit leaving angle is 
next considered in the same Figures. Comparing the 
two curves V = 0 and N = 0.2 (both for 5) per cent 
reaction, y = 1, = 0) we find that the influence of 
N is quite significant. For @ < 0.88, N reduces the 
exit distortion. The reverse is true ford >.88. Also, 
N reduces considerably the absolute magnitude of phase 
shift for all values of o. 

Since both NV and J/ exert strong influence on | 
as well as ¢,, the accurate prediction of the flow down- 
stream of an arbitrary inlet distortion consisting of 
many waves may be quite difficult even on the basis of 
this simple theory. 


PURELY CIRCUMFERENTIAL DISTORTION 
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Fic. 11. Influence of variation of loss and exit flow angle on 
ratio of total-pressure perturbation for purely circumferential 
distortion, 50 per cent reaction, y = 1, N = 0. 
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Fic. 12. Influence of variation of loss and exit flow angle on 


hase shift across rotor for purely circumferential distortion, 50 
Pp purely e 
per cent reaction, y = 1, N = 0. 


Although Eq. (36) will not apply for A, < up, it is 
interesting to examine the limiting case as \, — 0 and 
Ree Then Eq. (36) yields 

AU, (1 + a,b, +- Mk) (bo Nby,) + Mb,(Nk bs) 
AU (1 + — Nb:) — Marb, 


(42) 


which is areal number. Thus the phase shift along the 
circumferential direction decreases to zero as the dis- 
tortion velocity gradient becomes more severe span- 
wise’ than circumferentially. This must be true since 
for, the limiting case of purely spanwise distortion, 
the phase shift has been shown to be zero. [See case 
(1).] In particular, for 1J = 0 Eq. (42) gives 
AU2/AU, = (1 + aibi)/(1 + (43) 

and 
= (1 + aibs)(1 + + + (44) 
Case (4)—Combined Spanwise and Circumferential 
Distortion With un 

When the circumferential wave number \, is much 
less than the spanwise wave number y,,, corresponding 
to a distortion velocity gradient which is comparatively 
mild circumferentially and comparatively severe span- 
wise, w ~ u>>v. Here Eq. (13) or its equivalent (13’) 
must be fulfilled. We now have the four equations, 
(12), (13’), (19), and (28), for the four unknowns, 1, 
uo, AUs, and AV2. The solution is 


AU: 2(1 + abs) — Mb — 2b — 


= (45) 


For M = 0, this equation becomes identical to Eq. 
(43). Although both Eq. (45) and Eq. (42) reduce to 
Eq. (43) for M = 0, the value of u, and v, are not the 
same in the two cases. For the case represented by Eq. 
(45), = = (AU, — AU,)/2; these induced 
and distortion velocities in turn determine the exit 
flow angle which is not otherwise restricted. On the 
other hand, for the case represented by Eq. (42), u and v 
must vary in such a way that the specified relation Eq. 
(16) between the exit angle and the inlet angle be 
obeyed. Since our practical interest is in the ratio 
AU2/AU), the fact that u, v differs in these two cases 
is probably of secondary importance. 

We may now conclude that the solution represented 
by Eqs. (36)—(40) appears to be generally applicable for 
all distortion wave patterns except a purely spanwise 
distortion. Its applicability to a purely circumfer- 
ential distortion has been demonstrated by letting y,, = 
QO and ym,» = 1; Eq. (36) then reduces to Eq. (30). On 
the other hand, by letting A,, = 0 and »,,, = 0, it has 
been shown that Eq. (36) becomes identical to Eq. (45) 
at least for 14 = 0. Thus, although Eq. (36) is strictly 
true only for u», <A, it appears also to be a reasonable 
approximation for A, < um, thereby covering all dis- 
tortion wave patterns (except a purely spanwise one). 


Self-Induced Distortion 


In the foregoing we have evaluated the exit velocity 
distortion for inlet distortions of various wave forms. 
It might appear at first glance that if the inlet distor- 
tion is zero, the exit distortion would have to be zero. 
This is indeed the case under normal conditions. Under 
certain special conditions, however, there is a possi- 
bility of an exit velocity distortion even if there is no 
inlet distortion. In fact, there can be two different 
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Fic. 14. Self-induced combined spanwise and circumferential 
distortion, 50 per cent reaction, y = 1.0, N = 0, L,/L.z = 3r. 


types of such “‘self-induced”’ exit distortions—a purely 
spanwise type and a circumferential plus spanwise 
type. The latter also includes the purely circumfer- 
ential distortion as a special case, and, as will be 
shown, must travel or “‘rotate’’ with respect to the 
blades. 


Purely Spanwise Self-Induced Distortion 


This type of self-induced distortion can theoretically 
happen after either a rotor or a stator. It occurs when 
the determinant formed by the coefficients in the left- 
hand side of the pair of Eqs. (20) vanishes—i.e., 


M = [2 + (46) 


Since M varies from approximately zero to large values 
when the inlet angle 6; or its tangent 5; varies, there 
always exits a situation for which the above equation is 
satisfied for any given cascade. Fig. 13 shows the 
necessary value of / for such self-induced distortion to 
occur for 50 per cent reaction blading with y = 0 and 
y = 1 and with N = 0 (solid curves) and N = 0.5 
(dotted curves). However, as will be seen later, the 
onset of circumferential self-induced distortion will 
probably have the effect of suppressing this kind of 
spanwise distortion. 


Traveling Self-Induced Distortion; Rotating Stall 


Consider next the possibility of having self-induced 
distortions along both the circumferential and spanwise 
directions. Eq. (36) shows that for such a distortion 
to occur, we must have C= D=0); ie., 


a + (be Mazb, + 

— RN — M) 
a Asbo kb,)N M (a2 b;) 

1 + k? + a? 
This pair of equations can be used to solve for M and N 


in terms of the flow angles. It is clear that, except by 
sheer coincidence, M and N cannot both satisfy the re- 


0 (47) 


0 


Il 


1959 


quirements at the same flow angles. (The presence of 
Ymn does not help the situation since it can only assume 
certain discrete values.) Hence, a stationary self- 
induced distortion along the circumferential direction 
cannot happen. 

The possibility arises, however, of an unsteady, self- 
induced distortion. When there is no distortion in the 
inlet absolute flow, kU loses its significance as the pe- 
ripheral speed of the blades. Instead, kU’ is still the 
relative velocity between the (exit) distortion and the 
blades, but now becomes a dependent variable of the 
problem. From the relation a = 6b — k (which now 
becomes the definition for a), a is now the tangent of the 
main flow angle when viewed on a frame of reference 
which travels with the distortion. Eqs. (47) can now 
be used to solve k and M (or VN). The results are 


_ + Gabe) — Nbi) + — N) 
+ 
and is the root of the following equation: 
(1 + a2?) [Ymn?(1 + be?) + (mn? + — + 
+ 22”) + Ymn?)(1 +a2b2) = (49) 
This equation can be reduced to a cubic equation for k: 


Ek® + Fk? + Gk + H = 0 (50) 


M (48) 


where 


E = —2b, 

F = + + bo? + vmn?2(1 + N) 

G = + + + — + 
H = (1 + dy”) (bo? + Ndi? + ymn2(1 + N)] 


For the special case of a purely circumferential self- 
induced distortion, v,, = 1, and Eq. (49) yields an 
extremely simple solution for k: 


k = [(1 + + (1 + = 


(1 + N*)/2b; cos? B. (51) 


Substituting into Eq. (48), we obtain, after much 
simplification, 


M = [(1 + + — dibe) (52) 


Eqs. (51) and (52) have been obtained by the author 
in a previous report;® they are, respectively, Eq. (27) 
and Eq. (30) of reference 9. In fact, one can introduce 
a phase lag 6 in the response of the exit flow angle and 
loss coefficient to a variation of inlet flow angle by 
using, instead of M and N, Me’ and Ne.” The result 
as given in reference 9 is 


Me — b, tan 6 + N* sec 6 


2 cos? B2(b; + tan 4) 
M = [(1 + &*)/(1 + sec + 
[2k/(1 + bs*)]—N*(sin 6 + bo cos 6)} (54) 


These equations reduce to Eqs. (51) and (52) for 6 = 0. 

The necessary / as given by Eqs. (48) and (50) for 
self-induced combined circumferential and spanwise 
distortions is plotted in Fig. 14 for 50 per cent reaction, 
y = I, and VN = 0. 
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entire range of ¢ the necessary M/ for the purely circum- 
ferential type is lower than that for other types involv- 
ing both spanwise and circumferential self-induced 
distortion. This shows that the purely circumferential 
type will occur first and, hence, is the only type that re- 
quires further investigation. This simplifies matters 
enormously. For the purely circumferential type, 
calculations are made on 50 per cent reaction (Fig. 15) 
and 75 per cent reaction (Fig. 16) with y = Oand y = 1 
and with N = 0 (solid curves) and N = 0.5 (dotted 
curves). These results show that, for 50 per cent reac- 
tion with V = 0, the minimum values of M are around 
2at@ ~ 0.5 for y = 0 and around 0.66 at @ =~ 0.25 at 
¥ = 1. At other values of ¢, the necessary MM would 
be still higher. Now, even the smaller value M ~ 0.66 
(corresponding to an increase of @ of 0.023 per degree 
increase of inlet flow angle) is large enough to indicate 
that the cascade is near stall. For this reason the self- 
induced circumferential distortion can be identified as an 
incipient ‘rotating sfa//.’’ However, for the somewhat 
extreme case of N = 0.5, corresponding to very low 
solidity, the necessary / actually vanishes at ¢ = 0.14 
for Y = 1. The self-induced distortion at that region 
is then definitely not a stall phenomenon. It so hap- 
pens that @ = 0.14 is much too low for ordinary com- 
pressor cascades, and this may be the reason that the 
experimental evidence of such a non-stall type of 
traveling distortion has not been reported in literature. 

Comparing Eqs. (46) and (52), it can be seen that 
unless NV > 1, which is highly improbable, the necessary 
M for purely circumferential distortion is less than that 
for purely spanwise distortion. One may thus conclude 
that the purely circumferential type of self-induced dis- 
tortion always occurs first. The consequent nonuni- 
formity in the circumferential direction would in all 
likelihood tend to suppress the self-induced spanwise 
distortion. 
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Fic. 15. Self-induced circumferential distortion, 50 per cent 
reaction. 
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Self-induced circumferential distortion, 75 per cent 
reaction. 


Fic. 16. 


Velocity of Rotating Stall 


As mentioned before, kU is the velocity of travel of 
the self-induced exit distortion relative to the blade, 
irrespective of whether the blade is a rotor or stator. 
For the general case of distortion waves along both the 
circumferential and the spanwise directions, k can be 
evaluated from Eq. (50). For the simple case of a 
purely circumferential distortion, k is given by Eq. 
(53) with an assumed time lag and by Eq. (51) without 
a time lag. Now the purely circumferential distortion, 
when occurring at the onset of stall, has been identified 
as a ‘rotating stall’ and, as such, has been investigated 
by a number of researchers, both theoretically and 
experimentally. It is pertinent to compare the results 
in the present analysis with those of the other authors. 
This is done in the Appendix. 


Delay of Rotating Stall by Cascade Geometry 


Since the self-induced circumferential distortion or 
“rotating stall’’ is considered harmful to blades, it is of 
interest to find out if some cascade geometry might not 
be able to delay rotating stall better than others. The 
answer to this problem lies in Eq. (52), which gives the 
values of \/ that are necessary for rotating stall to 
make its appearance. If the necessary .\/ is very large, 
say approaching ©, then rotating stall may never 
occur. Unfortunately, this would require extremely 
low or high flow coefficients, as shown in Figs. 15 and 
16. For the practical range of compressor cascade 
geometry, rotating stall occurs at an inlet angle at 
which the actual \/ value of the cascade reaches the 
necessary ./ of Eq. (52). The delay of rotating stall, 
therefore, requires the simultaneous consideration of 
two problems—how to increase the necessary ./ and 
how to decrease the actual \/. In the following an 
attempt will be made to answer the first question. 
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DELAY OF ROTATING STALL 


LINES ARE "LINES OF UNDESIRABLE ¢ P 
CIRCLES ARE “CIRCLES OF NEUTRAL SOLIDITY 


FLOW COEF. > 
17. Delay of self-induced circumferential distortion or 
incipient rotating stall. 


(The second question can only be answered by com- 
paring the actual \/ and N values of different cascades, 
based on the results of cascade testing.) First con- 
sider the influence of V on the necessary V/ for rotating 
stall. (It is probably safe to assume that, other things 
being equal, lower solidity corresponds to higher NV and 
vice versa.) Now Eq. (52) shows that if bb. > 1, 
higher N will decrease the necessary M and, hence, 
promote rotating stall. The reverse is true if bby < 1. 
At b\b. = 1, the value of NV will have no influence. Ex- 
pressing in terms of the conventional design parameters 
¢, ¥, and R, it can be shown that the condition b,b, = 1 
is equivalent to 


+ (¥/4)? = 


If @ and y/4 are used to represent the coordinate axes, 
the above equation describes a circle with center at the 
origin and with radius equal to the reaction ratio R. 
(See Fig. 17.) These circles might be called ‘‘circles 
of neutral solidity,”’ in the sense that, if the combination 
of ¢ and y lie outside of the particular circle correspond- 
ing to the value of R of the cascade stage, rotating stall 
can be delayed by increasing N or decreasing solidity. 
On the other hand, if the design parameters lie inside 
of this circle then it is better to have a smaller NV by 
using a higher solidity. For 50 per cent reaction 
blading, the usual values of @ and y appear to lie out- 
side of the “circle of neutral solidity,’ hence lower 
solidity will most likely serve to prohibit rotating stall. 
On the other hand, for 75 per cent reaction, it is likely 
that some combinations of design parameters may lie 
inside the circle. For such cases, lower solidity might 
make the matter a little worse, somewhat contrary to 
usual expectations. 
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With the influence of solidity out of the way, we now 
may consider the value of VV with N = 0. Eq. (52) 
then gives 


M = (1 + = [62 + (R ¥/4)°1/4(R + ¥/4) 


With fixed R and y, we can differentiate ./ with respect 
to @ and set the derivative equal to zero. This will 
give us the worst possible choice of @ (corresponding to 
minimum necessary ./) for the given combination of R 
and y. It turns out that such undesirable values of ¢ 
can be described by the simple relation @ + (¥/4) = R. 
These are straight lines in the ¢@ versus y/4 plane, with 
intercepts equal to R at both coordinate axes. These 
lines might be called “lines of undesirable ¢” in the 
sense that, if the parameters ¢, y, and R lie on these 
lines, the resistance to rotating stall is the least. The 
obvious alternative is to stay as far away from these 
lines as possible in either side of them. Fig. 17 shows, 
perhaps unnecessarily, such a plot. It must be em- 
phasized that the values of ¢, y, and R are the ones 
corresponding to the operating (off-design) condition at 
which rotating stall is likely to occur; they may differ 
considerably from the design values. 


Appendix 


For the purely circumferential type of self-induced 
distortion, its velocity of propagation has been derived 
and given as Eq. (53) with an assumed time lag and 
Eq. (51) without a time lag. These results are com- 
pared below with some of the previous ones obtained 
by other investigators. 


Comparison With Previous Theories 
For the special case of N* = 0), Eq. (53) becomes 


| tan 6 tan 
2 cos? Bo(tan 8, + tan 4) 


This can be compared with Sear’s channel theory® which 
gives 


1 ees tan 6 tan @ 
cos” B[2 tan 8 + (1 — tan? B)tan 6j 


(In Sear’s theory, 6 implies 8;; however, no distinction 
is made between 6; and @.) On the other hand, Eq. 
(53) yields, for 6 = 0, the previous Eq. (51), 


k = (1 + N*)/2 cos? B, tan B, (51) 
This can be compared with Marble’s theory’ which gives 


(608° Be/cos? i) + N™ 

2 cos? tan 

If both N* and 6 are zero, the comparison is as follows: 
Present work K = 1/2 cos? tan (51a) 


Sears® and Marble’ k = 1/2 cos? 8; tan ~; 


To compare the present results with Stenning’s 
theory,® we take Stenning’s result for incipient stall 
(long wave length) : 
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TABLE 1 


Stenning, Kriebel, and Montgomery 
Stator, Solidity = 1.0 


Bi Be k 2k cos? Bo 
63.8° 47 .5° 0.63 0.58 
65.5° 49° 0.695 0.60 
67.2° 5p” 0.74 0.61 
68.5° 52° 0.72 0.54 

Stator, Solidity = 0.5 
64.8° 56° 0.78 0.49 
68.5° 57° 0.87 0.52 
Rotor 
67° 43° 0.85 0.91 
48° 0.85 0.76 
Costilow and Huppert 
Rotor 
53°-55° 25° 0.48-0.52 0.79-0.85 
0.46-0.50 0.75-0.81 


55°-57° 26° 


k = 2(1 — C,)/sin 28; ~ 1/cos? B tan p; 


if C, is taken approximately to be 1 — (cos? 8;/cos? 62). 
This gives a value of k equal to twice the value given by 
Eq. (5la). However, since V* in Eq. (51) lies between 
(0.1 to 0.5, the value of k given by Eq. (51) does not 
differ quite so much from Stenning’s results. (In 
another derivation assuming downstream mixing, 
Stenning obtained a different result which reduces to 
Eq. (51a) exactly for long wave length.) 

The result of Rannie and Marble* agrees completely 
with Eq. (51) for the case of N = 0. When JN is non- 
zero, there appears a discrepancy between Eq. (51) 
and the corresponding equation of reference 3. 


Comparison With Experiments 


Although a number of experimental results on stall 
propagation have been reported, relatively few would 
apply to the present analysis, which is based on an 
isolated rotor or stator and also restricted to incipient, 
small-amplitude stall. 

The data given in Table | for incipient stall propaga- 
tion in a single rotor or stator are gathered (or esti- 
mated if not directly available) from two sources— 
Stenning, Kriebel, and Montgomery® and Castilow and 
Huppert.'° 

Comparing these experimental data with Eq. (51), 
it appears that they follow the prediction reasonably 
well. In fact, these data can all be explained without 
the introduction of a time lag, for they lie within the 
curves of Eq. (51) N = Vand 0.5. (See Fig. 18.) The 
only exception is the case of the rotor of Stenning, et al., 
in which the experimental value appears to correspond to 
N* = 1. Perhaps by coincidence, a search of the rotor 
performance (8 versus §;) of reference 8 turns out that 
N* is indeed fairly large (around 0.8). 

In the paper by Rannie and Marble,’ a comparison 
between theory and experiments (Fig. 3 of reference 3) 
shows that their experimental values of k are higher 
by 10 to 15 per cent than their theoretically predicted 
values. The equation used in their theory, Eq. (7.1) of 
reference 3, appears to be the same as Eq. (51) of this 
analysis but with V* = 0. If one takes N* to be 0.10 
to 0.15, which are reasonable values for the cascade 
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Fic. 18. Comparison of theoretical velocity of propagation of 
incipient stall with experiments. 


of their investigation, there would be almost complete 
agreement between theory and experiment. 
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Integration of Some Thermal Differential 
Equations 


B. Saelman 

2 Weight Engineer, Lockheed Aircraft Corporation, Burbank, 
alif. 

May 13, 1959 


SYMBOLS 
Ty = boundary-layer recovery temperature, °R. 
T = skin temperature, °R. 
Te = skin equilibrium temperature, °R. 
Ta = ambient temperature, °R. 
To = initial skin temperature, °R. 
r = recovery factor 
a@ = acceleration Mach Number/sec. 
a0 = initial acceleration Mach Number/sec. 
Mo = initial Mach Number 
t = time, sec. 
= resistance, sec. (°R.)/B.t.u. 
Cn = capacitance, B.t.u./°R. 
m = heating rate, °R./sec. 
€ = emissivity 
o@ = Stephan-Boltzman constant 
q = heat flow, B.t.u./sec. 
hp = heat-transfer coefficient, B.t.u./ft.2(hr.), °R. 


M** OF THE EQUATIONS describing thermal relationships are 
solved by the finite-difference methods. However, there 
are some important equations in this category for which closed 
solutions exist. Some of these are presented below. 

For accelerated flight, we have the following relationships: 


a =a ct (1) 
T, = Tall + r[act — (ct?/2) + (2) 

The differential equation of the thermal relationship is 
(T; — T)/R = Cn(dT/dt) (3) 


Combining Eqs. (1)-(3) and integrating with the initial condition, 
T = 0 when? = 0, 
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where k = 1/RC, (5) 
jl 2r(ao? — Moc) rMy? 6racc Gre? 2ray Mol 
lk k ks ks 
When speed is constant, = = 0 and 
T = Ta(1 + 7M,2)(1 — e*) = T,(1 — e“**) (7) of the heat input is equal to the amount of heat radiated, we have 
During acceleration, the recovery temperature, 7, can sometimes hA(T, — Te) = AceT (20) 
be expressed Expressing = [T, — (T; — (21) 
Tr = mt (see Fig. 1) t < ti (8) Expanding [T7,; — (7; — T;)]*, neglecting terms involving powers 
Substituting Eq. (8) into Eq. (3) and integrating of T, — T¢ greater than unity, and solving, 
T = (m/k)(kt — 1) + (m/k)e™ t<t (9) Te = [T, + + 40(€/hy)T;*] (22) 
and T; = (m/k)(kt; — 1) + (m/k)e* (10) This approximation is valid when 7, — T, is small. 
and T, = T; constant for t>h (11) 


Substituting Eq. (11) into Eq. (3) and integrating 
(m/k)A (12) 


The temperature difference can be expressed 


Ta = T, —T = (m/k)(1 — (13) 
This isa maximum at t = 4; 
Tamaz = (m/k)(1 — e~“) (14) 
Introducing the radiation term 
q = AoeT* (15) 
q = (T; — T)/R = C,(dT/dt) + AoeT (16) 


The integral of Eq. (16) is 


+ kT — kT, 
AdT BdT (CT + D)dT 
T-a e7?+/T +2 
where a and Bb are real roots of 
(Age/C,)T* + RT — kT; = 0 
t= — (A log(T — a) + Blog (T — b) + (C/2e) X 
log eT? + fT + g| + {21D — (Cf/20)|/W4eg — X 
tan! [(2eT + f)/W4eg — (18) 
Since this expression involves T as the independent variable, the 
radiation term is usually neglected and integration of the differ- 


ential equation leads to the convenient expression, given below, 
with the temperature expressed as a function of time: 


= T, —(T, — (19) 


Neglecting the capacitance term in Eq. (16)—i.e., assuming all 
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On an Exact Solution of Incompressible, 
inviscid, Hypersonic, Stagnation-Point Flow 
for a Sphere* 


Koon-Sang Wan 
Research Associate Professor of Applied Mechanics, Polytechnic 
Institute of Brooklyn, Freeport, N.Y. 


April 10, 1959 


’ | ‘HE PURPOSE OF THIS NOTE is to present an exact solution for 
an inviscid hypersonic flow at the stagnation point of a 
sphere. The density p of the flow in the ‘‘shock layer” is assumed 
to be constant and is equal to p/p, = 1/k (see references 1 
and 2) where the value of k depends on the free-stream Mach 
Number /,, and the ratio of the specific heats y = C,/C,. 
Under these assumptions, the governing equations can be writ- 
ten as, in spherical coordinates (see Fig. 1), 

(1/r?)(0/dr)(r2u) + (1/r sin 0)(0/06)(v sin 6) = 0 (1) 
u(Ou/dr) + (v/r)(Ou/00) — (v?/r) = —(1/p)(Op/dr) (2 
u(dv/dr) + (v/r)(dv/d0) + (uv/r) = —(1/rp)(dp/d0) (3) 

The boundary conditions are: (1) at r = a, 0 = 0°; u = 0, 
vy = 0; (2)atr =a,u = 0, and (3) atr = a + 6(0) =a + 
i.e., near the stagnation point, 
u = u.,{—[1 + (k — 1) sin? B] cos @ + 
(1 — k) sin B cos B sin 6} => —ku,, cos@ (4) 
v = (1 + (k — 1) sin? 8] sino + 
(1 — k) sin B cos B cos 6} ~u,, sin@ (5) 
b= Pao t+ — k) sin? 
Po — k)cos?@ (6) 


by using the oblique shock conditions and where 
= 90° — 6 + [1/(a + = 90° — 
The solution of this problem can be obtained by assuming 
u =f(r)cos0, v = g(r) sin 6 (7) 
Substituting these into Eq. (1) yields 
2g = —(1/r)(d/dr)(r°f) (8) 
By cross-differentiating Eqs. (2) and (3) to eliminate » and using 


* The results obtained here are identical with those obtained by Lighthill,* 
but the approach used here is slightly different. 
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relation (8), one can obtain another relation between f and g as 
r(d*g/dr?) + 2(df/dr) = 0 (9) 
The solutions of f and g are, therefore, G 
f(r) = A, + Aor? + (A;/r%), 
g(r) = — 2Azr? + (1/2)(A3/r?) 
u = [A; + Aor? + (A3/r)] cos 0 (10) 
v = [—A, — 2Aor? + (1/2)(A3/r*)] sin (11) 


Hence, 


and 
p/p = As + (1/2) [AiAor? + (A1A3/r3) + Aor? + 7A2A3(1/r) — 
(1/4)(A3?/r)] — (1/2)[BAiAor? + + 
2A2*r4 + 9A2A3(1/r) + (3/4)(A32/r*)] cos? @ (12) 
The unknown constants A;, A», A3, Ay, and 6* are now deter- 
mined by the boundary conditions and obtained as 
= A,/u,, =(1 — 4k/3k 
A,(a + 6*)?/u,, = —(1 — k)?/5k 
A; = A;/[u.(a + = 2/1 — k)(6k — 1)/15k (18) 
Ay = k(po/Po) — + + 
— 


i] 


1 
A, 


and 
—3(1 — k)?o® + (5 — 20k)o* + 2(1 — k)(6k —1) = 0 (14) 


an equation determining the value of o = a/(a + 6*) and, hence, 
5*. Since & is related to M,, therefore u, v, and p within the 
shock layer and near the stagnation point are completely deter- 
mined for a given M,, from Eqs. (10)—(14). 

It is of interest to compute the values of 6*/a vs. M,,. The 
results of such computations are shown in Figs. 2 and 3 together 
with results from experiments and theory by Li and Geiger for the 
purpose of comparison. They show that the present results check 
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gases). 


more closely with the experiment than those of Li and Geiger. 
This is expected since the effect of curvature is included in the 
present solution, as was also pointed out by Li and Geiger. The 
pressure coefficient at r = a has also been computed vs. 6. The 
results are in good agreement with those from experiments and 
theory of Li and Geiger near the stagnation point. 
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Comments on ‘‘The Heating of Slabs With 
Arbitrary Heat Inputs”’ 


William Squire* 
Southwest Research Institute, San Antonio, Texas 
May 19, 1959 


I A RECENT COMMUNICATION, Goodman! described the applica- 
tion of a ‘‘heat-balance integral’ method to heat conduction 
problems including those with nonlinear boundary conditions. 
The basic advantage of the method is that approximate solutions 
of a partial differential equation are obtained by solving an 
ordinary differential equation. This note shows that an ap- 
proximation method due to Shvets? gives analogous results and 
points out the need for further investigation of these methods. 
The basic problem to be considered is the partial differential 
equation 


07u/O0x? = Ou/Odt (1) 
subject to the initial condition 
u(x,0) =0 (2) 
and the boundary conditions 
u(0, t) = 0 (3a) 
and ou(0, t)/Ox = —f(z, t) (3b) 


where is the surface temperature, ¢). 
The solution is approximated by 
k=n 
w= > (4) 
k=0 
n being an integer > 1 and & an integer > 0. The partial solu- 
tions, u,, satisfy the equation 


* Formerly of Bell Aircraft Corporation, Buffalo, N.Y. 
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= Ouj,—1/Ot (5a) 
and the sequence is started by 
= 0 (5b) 


This gives rise to polynomials in x so u” cannot satisfy Eq. (3a). 
Therefore this boundary condition is replaced by 


u"(6,t) = 0 (6a) 
and the additional condition 
(Ou"/Ox)(6, t) = O (6b) 


is imposed to determine 6(/), the distance the heat penetrates the 
slab. Two boundary conditions are needed to determine the 
integration ‘‘constants”’ of Eq. (5) and two conditions are needed 
to determine 6 and z as functions of ¢. There are three possible 
choicest which are shown in Table 1 where the corresponding 
uo, U1, and equations for 6 and s are also listed. 

We cannot present a detailed discussion of the solutions of the 
equations for 6 and z here but we will show that marked differences 
exist in the predictions of the different approximations. For the 
case where f is a function of ¢t only, Goodman obtains 


= (12F/f)!/? = 32/f (7) 
where 
t 
F(t) = 
The supplementary equations of Case 2 can be integrated to give 
= (4F/f)!/2 (8) 


which is similar in form but with a less accurate numerical 
coefficient. On the other hand, Case 2 gives 


(2 = 2F/s (9) 


Eq. (9) appears to be superior from a logical point of view since it 
gives 6 as an expression which is inherently positive while Eq. (7) 
changes sign when f or z change sign. Dr. Goodman has foundt 
that his method is not suitable for pulse-type inputs and seems to 
break down when f’ changes sign. 

On the other hand, if we consider the case where 


fo) =t (10) 
t Shvets’ paper deals primarily with boundary-layer problems, and the 
treatment of heat conduction is very concise. The author interpreted it as 
Case 1. Dr. Goodman, to whom the author is indebted for an interesting 
correspondence, interpreted it as Case 2. This led the author to consider all 
three possibilities. 
t Private communication. 


TABLE 


EQUATIONS 
INTEGRATION SUPPLEMENTARY PARTIAL FOR z AND 8 
CASE CONDITIONS CONDITIONS { SOLUTIONS _ (FIRST ORDER) 


Tan 


au" 82! ea! 
Up (0,1) =O (8.10 2t=8'z + 82! 
2 
uy (8,1) =0 
é 
uy (8,1) - (x3 
M up l0,t)=z2(t) (8,1) 0 Uo z= 
n ¢! ¢! 
up (0,1)= 0 0? 


du, 
Fe (0,1) 


(0,1) 


k IS AN INTEGER 20 n IS AN INTEGER 21 


SIs 


note: f'= 


the exact solution is 


9 


= {n!/[n + (1/2)]}!} (11) 
Goodman's method gives 
[4/3(1 + + 4/2) (12) 
while Case 1 gives 
s = + + + 2) (13) 


A comparison shows Goodman's method is in closer agreement 
with the exact solution. 

It, therefore, would be desirable to investigate the various 
approximations further in order to determine their relative 
advantages for different types of problems and to put them on a 
sound mathematical basis. 
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Note on the Calculation of Equilibrium 
Skin Temperatures 


Hans C. Vetter 

Aerodynamics Engineer, Douglas Aircraft Company, !nc., 
El Segundo, Calif. 

May 21, 1959 


—— NOTE presents a method for the rapid estimation of 
equilibrium skin temperatures that result from aerodynamic 
heating in both the real-gas and ideal-gas regimes. Heretofore, 
the calculation of equilibrium temperatures has relied on a graphi- 
cal or analytical determination of the roots of a fourth-order 
equation. The present method reduces the problem to the read- 
ing of a design chart for ideal-gas flow and to a quickly converging 
iteration method using the same design chart for real-gas flows. 


IDEAL GAS 
The heat-balance equation for an ideal gas is 
h(T, — + aG = (1) 
The first term on the left of Eq. (1) is the aerodynamic heating 
input, the second is an irradiation term, and the surface-radiation 
term appears on theright. By defining 


= h/eo 
Eq. (1) may be rewritten, as follows: 
+ = (Tr + aG/h)/r (2) 


The roots of Eq. (2) may be calculated for all possible values of 
the right-hand term, and the results plotted to form a design 
curve for the determination of equilibrium temperatures. Actu- 
ally, it is much easier to plot the right-hand term versus 7,/r. 
A plot of Eq. (2) is shown in Fig. 1. 

A temperature-dependent heat-transfer coefficient may be 
handled by iteration, but it is often sufficient to make the initial 
evaluation of h at a temperature somewhat lower than the re- 
covery value. 

REAL Gas 

The appropriate heat-balance equation for aerodynamic heat- 

ing with real-gas effects but neglecting hot-gas radiation is 
Cup V(tr — te) + aG = eoT,! (3) 


Here, Cy is the Stanton Number, p the density, V the velocity, 
and i the enthalpy. Now, the following definitions are made: 


— 3 
lica- 
‘tion 
ions. 
ions 
an 
ap- 
and 
Ss. 
itial 
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(2) 
3a) 
3b) 
(4) 


h = Cup 


te = Cy Te + Ate 
T,’ = 1,/¢p, — Ate/Cp, 


where c,, is the ideal-gas value of the specific heat. Using these 
definitions, Eq. (3) may be written as follows: 


h'(T,’ — T.) + aG = eoT 4 (4) 


Eq. (4) has exactly the same form as Eq. (1), with the exception 
that 7,’ is a function of 7... Hence, Eq. (4) may be put into 
the generalized form of Eq. (2) using the primed quantities in 
place of the unprimed quantities. By assuming Az, to be zero, 
a first estimate of 7, may be calculated from the design chart 
that may then be used to calculate a second approximation to 
T,’. By using this value of 7;,’, a better estimate of 7, may be 
made, ete., until the desired accuracy is achieved. A temper- 
ature-dependent Stanton Number may, of course, be included 
in the iteration. 

If real-gas effects are not large (Ai, << i,), but not negligible, 
one iteration is sufficient to determine the equilibrium temper- 
ature. In addition, it is often sufficient to use the following 
approximation: 


T, = T,, + AT, 
hy’ = hy’ + Ah,’ 
AT, = + (OT./OT,’)AT,’ 


whence 
AT. & {[(ir/C,,) — — /(ho’ + 


where the subscript 1 refers to an evaluation at 7, = 7,,, and ho’ 
is the initial estimate of hk’. The above approximation is based 
on the observation that the equilibrium temperature, as calcu- 
lated from Eq. (4), differs from that calculated from Eq. (1) by 
only a small difference in the ‘‘effective’’ recovery temperature 
(AT,’ = —Ai./C,). 
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On the Integral of Moment of Pressure 
in Supersonic Flow* 


Lu Ting 

Research Associate Professor, Polytechnic Institute of Brooklyn, 
Freeport, N.Y. 

May 21, 1959 


N THIS PAPER, the integral of the moment of the disturbance 
pressure along the line of the interaction of a Mach plane 


* This research was supported by the USAF through the Air Force Office 
of Scientific Research, Air Research and Development Command, under 
Contract No. AF 49(638)-217, Project No. 9781. 

The author wishes to express his sincere thanks to Prof. Antonio Ferri 
for many stimulating discussions and suggestions. 
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with the cylindrical surface whose generator is parallel to the 
undisturbed stream, is shown to be related to the integral of 
the moment of the prescribed normal velocity along the same 
line. It is derived by the generalized integral relationship which 
relates the integral of disturbance pressure to that of the pre- 
scribed normal velocity.' 

As shown in Fig. 1, a cylindrical surface y = f(s), s = g(s), is 
placed in a supersonic stream with its generator parallel to the 
x-axis, which is the direction of the undisturbed supersonic flow 
with velocity U and Mach Number 1. q,(x, s) represents the 
prescribed small normal velocity on the cylindrical surface with 
gn = 0 forx <0. It will be assumed that for any given Mach 
plane, x + B(ycosw + zsinw) = X, there exist two finite values 
of Ki(X) and K.(X) with K; > Ko, such that ahead of the Mach 
plane, the disturbance potential due to g, is confined inside the 
region K,; >z> Ko; i.e., 


glx < X — Biycosw +z sinw), y, zs] = 0 (1) 


for > K,(X) and This assumption is 
usually valid for a finite range of w. The generalized integral 
relationship states: 


f fe y, sn-@dS = (ov /B) ff a dS (2) 


where p(x, y, 2) is the disturbance pressure; n is the unit vector 
normal to the cylindrical surface, ® = cos wj + sinwk; j, k are 
the unit vectors parallel to y and s axis; and S represents the 
domain of integration on the cylindrical surface ahead of the 


Mach plane. Eq. (2) can be written as 


ds dx p[x, f(s), g(s)](my cos + nz sin w) = 
x*(s 


where x*(s) = X — B[f(s) cosw + g(s) sin w]; 5; and sz are the 
are lengths on the cylindrical surface corresponding to the points 
s = K, and z = Kz, respectively; and m,, mz are the y and z 
components of the normal n which is a function of s only. Dif- 
ferentiation of Eq. (3) with respect to w yields: 
Se 
Bf ds p(x*, f, g)(m-@)(f sin w — g cos w) — 
Si 
ff p(x, f, g)\(n, sin w — n, cos w)dS = 
eu ds q,(x*, s)(f sin w — g cos w) 
$1 
which is equivalent to 
ds(n-@)(y sin w —zcosw)p = (1/B) X 
Sf, (n X @-i)pdS + (ou/B) ds (ysinw — 2cosw)gn (4) 


where the line integrals are carried out along the line 7, which 
is the line of intersection of the Mach plane (XY, w) and the cylin- 
drical surface. 
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Fic. 2. Planar wing at yaw with a straight trailing edge. 


Eq. (4) is the general relationship which involves the integral 
of moment of pressure. Eq. (4) isnow applied toa planar system. 
Fig. 2 shows a wing planform with all edges supersonic and a 
straight trailing edge. The wing is yawed at an angle A. The 
total lift and center of pressure can be calculated readily by virtue 
of Eqs. (2) and (4). For the planar system withn = jand y = 0, 
Eqs. (2) and (4) become 


cose ff pas = (5) 
and 


Beoste fi, spdz = -sine [feds + ov fi sand: = 
7 
pl [ + fi. | (6) 


For w = —arc sin (tan A/B), the line T is parallel to the trail- 
ing edge and S is the area of the upper surface (or lower surface ) 
of the wing ahead of the line 7. For X = Xo, the line 7 coin- 
cides with the trailing edge, and the total lift acting on the upper 
surface of the wing is 


L+ = p* dA = —(pU/B cos w) Sf, dA (7) 


where the superscript + refers to the upper surface. The loca- 
tion (z, x) of the center of pressure acting on the line T is ob- 


tained from Eq. (6): 
+o lo bt + de — tz 
2X) = (2 is —tanw X 


and x(X) = X + ztan A where the relationship 


cose p* dz = (pU/B) f, Qn dz 


is deduced from Eq. (2).! The location of the center of pressure 
on the upper surface of the wing is simply 


x,t = | — x(X) p*) dX/L* 
as) ax / ff, ant dA 
2(X) ax/ ff, (9) 


Similar results can be obtained for the lower surface. Additional 
applications to planar and nonplanar systems are presented in 
reference 3. 
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Column End-Fixity Versus End-Restraint 


Paul E. Sandorff 
Associate Professor of Aeronautics and Astronautics, 
Massachusetts Institute of Technology, Cambridge, Mass. 


May 21, 1959 


¥ A RECENT REPORT published here, Neffson' has made a 
creditable presentation of a problem that recurs repeatedly in 
practice and in the literature. The earliest solution of this prob- 
lem with which the writer is familiar is that of Osgood,? who, em- 
ploying the same fundamental beam-column theory, derived a 
method and a chart whereby the critical load could be expressed 
from the end-moments applying to an elastically-restrained test 
column. Osgood references earlier work by Prager and Zimmer- 
man, which was limited to consideration of elastic columns. 
Wozniak* presented a solution to the inverse problem of deter- 
mining the end-moments, in the presence of a lateral load, when 
the end-fixity coefficient is established. Most certainly unknown 
to Neffson is reference 4, which employed in part the identical 
derivation as reference 1 and presented the suggested curves of 
end-fixity coefficient vs. end-restraint factors, where the latter 
are established by the relative rigidities of the column and the 
adjacent restraining structure. Deserving of mention is the un- 
published report of Ensrud, Kastan, and Wignot® on working 
methods for all aspects of beam-column structures, including 
those continuous over many supports and with various types of 
lateral loads and restraint conditions. There are undoubtedly 
others whom the writer regrets being unable to credit. 

Hitherto unpublished is the much improved and extended 
chart expressing column end-fixity in terms of end-restraints, due 
to C. H. Heilbron, presented herewith as Fig. 1. 

A typical problem which is greatly simplified by use of such a 
chart is the determination of a proper test length, at, say, c = 3.5 
for a structure which, in service, is stated to have ‘‘end-fixity con- 
ditions corresponding to c = 2.0 at one end and toc = 1.5 at the 


other.”” Such improperly stated specifications may be trans- 
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lated, ‘‘if end-restraints were symmetrical and equal to a, c 
would be 1.0; if both were equal to a2, c would be 2.0... From 
Fig. 1, equal end-restraint factors corresponding to c = 2.0 are 
@L/EI = oL/EI = 3.4; corresponding toc = 1.0, aL/EI = 
= 0. The actual end-fixity which corresponds to a,L/EI 
= Oand amL/EI = 3.4, isc = 1.44. The proper test length is 
therefore, /3.5/1.44 times as large as the actual structure. 

The chart extends into the regions of negative end-restraint, 
so that it may be used to advantage when destabilizing influences 
are exerted by the adjacent structure. It is, therefore, suitable 
for the solution of multispan columns as suggested in reference 4. 
A typical problem concerns a column of constant section continu- 
ous over a central support, but, for example, with one bay half 
the length of the other, and where all supports are assumed to 
furnish restraint corresponding to c = 1.5. Simultaneous insta- 
bility requires that cap/cne = (L1/L2)? = 1/4. By trial and error, 
a solution is quickly found whereby the proper increment of 


stiffness is exchanged at the central support; it is agL/EJ = 1.4; 
a-pL/EI = 1.4 — 3.66 = —2.2, a+,L/EI = 1.4 + 2(3.6) = 8.6, 


a-L/EI = 1.4, and cap = 0.50+, cre = 2.01. 
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Comments on ‘‘A New Solution to the Laminar 
Boundary-Layer Equations”’ 


M. Z. v. Krzywoblocki 
Professor, University of Illinois, Urbana, Ill. 
May 21, 1959 


Hs proposed a generalization of Gértler’s solution of 
plane and steady laminar boundary-layer equations for 
a certain class of outer pressure distributions. The remarks 
below refer to both these propositions, and, as a matter of fact, 
to many similar and analogous propositions of this kind. The 
equations of a boundary layer are nonlinear partial differential 
equations which are usually associated with two-curve boundary 
conditions. In the theory of differential equations, one dis- 
tinguishes two main processes: (1) formal construction of a 
solution satisfying the prescribed boundary conditions; (2) proof 
of existence and, possibly, of uniqueness of a solution—i.e., proof 
that a solution actually exists and satisfies boundary conditions, 
regardless of its form (closed or derivable by means of a limiting 
process). In both processes, one may attempt to derive a solu- 
tion and a proof in an infinitesimal neighborhood of the boundary 
curve (local solution or solution in the small) or in a finite region 
(solution in the large). Up to the present time, there do not 
exist any proofs of the nature of process (2) referring to the 
boundary-layer equations in question. Consequently, all the 
attempts and efforts referring to these equations can only be 
classified as processes belonging to category (1). Since there is 
no existence proof, one cannot investigate whether a solution 
constructed by some method is exact, convergent, approximate, 
or asymptotic (this may be even divergent). A possible com- 
parison may be made with physical (experimental) test results. 
The classical mathematical approach to the field of differential 
equations is one of demonstrating that the problem in question 
has at least one solution (sometimes proving uniqueness). In 
many cases, an investigator may well seek a series or numerical 
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solution even though he is not certain it is obtainable in the 
manner chosen (or at all). The purpose may be an application 
or it may be one of experimental arithmetic—i.e., a calculation 
which may help in an effort to guess a correct statement of the 
proposition in question. If one believes that there is a prob. 
ability of a success and if one can verify that his final result is an 
adequate approximation to a solution, then there may be no rea- 
son to insist on furnishing a firm proof. But the status of the 
entire question of the so-called new and old solutions, supposedly 
known as exact or rigorous, etc., should be always explained and 
discussed, and it is believed that these few remarks may help to 
understand it a little better. 

An interesting contribution to the problem in question in the 
incompressible-fluid flow domain was furnished by Rheinboldt 
who proved that if in a certain region [0 S x S$ m;0 Sy S @] 
the stream function y is an analytic function of x, then the 
boundary condition at infinity is always automatically satisfied 
and does not need to be taken into account. But to this end the 
function y must be representable by an absolutely convergent 
power series in x. In the existing texts and collective works 
(Princeton Series, Oxford Series, etc.) on boundary-layer theory, 
the description and explanation of these fundamentals is usually 
neglected. 
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Further Remarks on Matrix Interpolation of 
Flexibility Influence Coefficients} 


William P. Rodden* 
Research Specialist, Vibration, Flutter and Acoustics, North 
American Aviation, Inc., Los Angeles, Calif. 


May 21, 1959 


I CARRYING OUT an aeroelastic analysis by collocation methods, 

it is assumed that both structural and aerodynamic influence 
coefficients are known. The aerodynamic influence coefficients 
have been defined! as relating the forces acting at the control 
points to the deflections of the same control points. Since the 
structural influence coefficients normally are determined for 
loads acting at the panel points of the structure, it becomes 
necessary to unify the control points of both the structural and 
the aerodynamic systems in order to obtain a consistent aero- 
elastic analysis. It is the nature of aerodynamic theories that a 
transformation of control points is not convenient, whereas the 
transformation of the structural control points can be carried out 
readily. 

We will generalize the results of Schmitt? in order to eliminate 
the dependence on the least-squares interpolation method. We 
denote the structural net by p and the aerodynamic net by R. 
We are given the structural matrix [a,] such that 


{hp} = {Fp} (1) 
and we seek an equivalent set such that 

{he} = [ax] | Fe} (2) 
The problem is to relate { F,} to {F,}. We choose to make the 


relationship structurally equivalent rather than statically 


equivalent because of the greater accuracy of the outboard 


t The present analysis is an outgrowth of a private communication with 
Dr. Alfred F. Schmitt, and the author wishes to acknowledge his very 
helpful and stimulating comments. 

* Now Engineering Specialist, Dynamics Branch, Norair Division of the 
Northrop Corporation, Hawthorne, Calif. 
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deflections which are more important in an aeroelastic analysis. 
By structurally equivalent we mean that the set of forces | F,} 
produces the same deflection mode {hp} as does the original set 
\F,}. Therefore, there is some matrix [C,.] that defines the 
structural equivalence. 


= [Cox] Fe} (3) 
By Maxwell’s Law of Reciprocity 
the} = [Cup] | Fo} (4) 
where [C;,] is the transpose of [C,x.]. From Eqs. (1) and (4), 
= [Crp] hp} (5) 
and from Eqs. (3) and (5), 
{he} = [Crp] Fe} (6) 
By identifying Eqs. (2) and (6), we find 
[ae] = [Ckp]lap]* [Cpe] (7) 


To determine [C,x], we assume a deflection interpolation 
matrix can be found. 


{he} = {hp} (8) 
By identifying Eqs. (5) and (8), we find 

Zip] = [Ckp] [ap]~ (9) 
from which [Cpx] = [ap] [Zp«] (10) 


Substituting Eq. (10) into Eq. (7) vields the desired form of the 
transformation. 


fax] = [ap] (11) 


The result Eq. (11) is that obtained by Schmitt if the interpola- 
tion matrix of Eq. (8) is found from the least-squares method. 
However, Eq. (8) shows that there are a number of methods obtain- 
ing the elements of [J;,]. These methods can be divided into two 
classes: interpolation-in-the-large in which the entire surface is 
fitted with a single set of functions—e.g., Schmitt suggests that 
product combinations of the normal modes of a uniform canti- 
lever beam provide satisfactory functions for plate-like wing 
surfaces, and interpolation-in-the-small in which each small 
region is fitted with functions appropriate to the significant 
boundary conditions and structural discontinuities in the region 
(in this manner, the elements of [J;p] are found a row at a 
time). It appears that the latter choice is the more reliable for, 
although a practical surface is uniform for the most part, there 
are always a sufficient number of structural discontinuities to 
be bothersome. The interpolation-in-the-small can be carried 
out using double polynomial functions (x"y") either by the 
Lagrangian method (which has been found to be grossly in- 
accurate for interpolation-in-the-large) or by the least-squares 
method. The least-squares method is more general since it 
includes the Lagrangian fit as a limiting case. 
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Visual Supersonic Flow Patterns by Means of 
Smoke Lines 


V. P. Goddard, J. A. McLaughlin, and F. N. M. Brown 
University of Notre Dame, Notre Dame, Ind. 
May 26, 1959 


TT CONSIDERABLE SUCCESS enjoyed here at Notre Dame in 
extending the useful top speed of three-dimensional smoke- 
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tunnels! from a bare 20 ft. per sec. to more than 200 ft. per sec. 
resulted in a doubt that any upper limit of speed for smoke 
tunnels existed. The radical increase in useful speed in the sub- 
sonic tunnel was achieved rather simply by severe reduction 
ratios, heavy antiturbulence screening, and introduction of the 
smoke upstream of this screening. The same technique of 
design was applied to the supersonic tunnel, and the doubt was 
confirmed; smoke pictures can be taken at least as high as 
1,325 ft. per sec. 

The techniques for producing and photographing supersonic 
flow patterns have not yet been fully perfected. However, they 
are at a point which enables a realization of the potentialities of 
this type of aerodynamic study in the supersonic flow regime. 

It is the purpose of this article to present and describe these 
techniques. Pictures are included of flow patterns about a 
wedge and a sphere (see Figs. 1 and 2). 


TUNNEL 

The Notre Dame supersonic pilot tunnel was designed by the 
Foelsch method.? It is an indraft type tunnel capable of con- 
tinuous operation at Mach 1.38. The working section is ap- 
proximately 4 in. long with a cross-sectional area of 6.25 sq. 
in. The entire side walls are of 3/4 in. lucite to enable 
visualization and schlieren setups along the full length of the 
nozzle. 

The reduction cone was designed by the Smith and Wang 
curves.* It has a reduction ratio to the throat of the tunnel of 
96 to 1 and is 3 ft. in length. There are 6 antiturbulence screens 
of fine-mesh nylon at the upstream end of the reduction cone 
spaced 3/4 in. apart. 


Fic. 1. Smoke-flow pattern about a wedge with shock 
shadowgraph superimposed. Mach Number = 1.38. 


Fic. 2. Smoke-flow pattern about a sphere with shock 
shadowgraph superimposed. Mach Number = 1.38. 
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SMOKE SYSTEM 
The smoke is produced according to techniques devised by 
Brown.‘ Wheat straw is coked, and the resulting smoke is cooled 
and filtered before introduction. 


LIGHTING AND PHOTOGRAPHY 

The light source for both shock shadow and smoke filaments was 
a short-interval flash lamp synchronized with the shutter of a 
modified schlieren camera. The effective exposure was about 20 
microsec. 

To make pictures of the shock shadow, the short-interval lamp 
was substituted for the usual steady light source. A cover with 
a pin-hole opening was placed over the light to produce a point 
source. The knife edges were omitted, and the image of the 
point source was focused on the shutter diaphragm. 

The pictures of the smoke lines were made with the same 
modified schlieren camera, but two short-interval lamps were 
placed behind and upstream of the working section. Both lamps 
were aimed at the model. The angle of the beam with respect 
to the back wall of lucite was just larger than the critical angle of 
that material. Tungsten 160 film was used, and the desired 
contrast was obtained by underexposing and overdeveloping. 

The accompanying pictures were made through the two 
negatives—the smoke-filament negative and the shock-shadow 
negative. This is justifiable because the negatives were identical 
in scale, a steady-state flow condition existed, and both negatives 
were made with smoke flowing. The random marks on the 
pictures are due to scratches on the lucite tunnel walls. 


FUTURE PLANS 


The current pictures are, perhaps, acceptable, but it is believed 
that they can be considerably improved. To this end, a new 
nozzle is under construction which permits introduction of the 
light through a lucite plate sandwiched in the lower nozzle 
block. This will permit full utilization of the polarization prop- 
erty and, at the same time, diminish the illumination of the 
background. 

Since the level of turbulence appears to be a function of the 
reduction ratio and the screening, the optimum will be sought. 
The reduction ratio and the number of screens and screen 
spacing will be varied independently to this end. It is hoped an 
optimum combination can be found to enable taking pictures of 
the quality now being obtained with the subsonic smoke tunnel 
here at Notre Dame. 
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The Compressible Laminar Boundary Layer 
With Constant Injected Mass Flux at the 
Surface 


Alan Q. Eschenroeder 
Cornell University, Ithaca, N.Y. 
May 26, 1959 


i of the steady compressible laminar boundary layer 
on a flat plate with fluid injection have been done by Low! 
for injection velocity proportional to the reciprocal of the square 
root of the distance from the leading edge and by Lew and 
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Fanucci? for uniform injection velocity. Certain problems jp 
mass-transfer cooling with combustion require a compressible 
solution for uniform mass flux from the surface because the wide 
variation of fluid density precludes the use of constant injection 
velocity as a boundary condition when fluid is supplied uni- 
formly to a porous wall boundary. 

Under the assumption that the product of density and vis. 
cosity is a constant, the uniform mass flux case can be solved in 
terms of an existing incompressible solution by Lew and Fanucci.! 
Basically, their transformation procedure is followed with re. 
definitions of variables. The boundary-layer equation for the x- 
component of momentum after the von Mises transformation 
reduces to 


ou'/dx’ = (0/dy’)[u'(du'/dy’)] (1) 


u’ is the nondimensional downstream velocity (referred to the 
free-stream velocity U,), x’ is the nondimensional downstream 
distance from the leading edge (referred to a typical downstream 
length L), and y’ is a nondimensional stream function defined by 


p’u’ = /dy’); = —(1/*~/ (2) 


where p’ is the nondimensional density (referred to the density 
at free-stream conditions p;), y’ is the nondimensional normal 
distance (referred to L), v’ is the nondimensional normal velocity 
(referred to U;), and Re, is a Reynolds Number p,U,L/u. The 
boundary conditions become 


= 1; a(x’, = 1; = 0 (3) 


where f’(x’) is given by 


= —V Re, 


The subscript 0 refers to conditions at the wall surface. Further 
transformations following the method of reference 3 with con- 
sistent redefinitions of variables lead to the equation 


V + — (WV'/r') + - 
= (85) 


, 


dx’ (4) 


where V’ = 4u’2 and the two independent variables are related to 
the original variables by the equations 


= 2[ —(povo/ pr (6) 


r’ — y 
and 20’ = —(povo/m) (p/pi)dy (7) 


The boundary conditions governing Eq. (5) are 
V‘(o’, 0) = Oand Vo’, ~) = 4 (8) 
Eq. (5) with these boundary conditions is formally identical to 
the problem solved numerically in reference 3; therefore, the 
tabulated results of that reference can be utilized for the com- 
pressible case with uniform mass flux if the definitions expressed 
in Eqs. (6) and (7) are used. As a check, it can be shown that 
the normal distance in reference 3 is given by 


ff.” 


according to the familiar Dorodnitzyn-Howarth transformation, 
which should be implicit in the present method. 
The friction coefficient c; defined by 


Cr = wo Ou/dy))o/(1/2) ai Ui? (9) 
can also be found. For various values of o’, Lew and Fanucci’ 
have tabulated a starting coefficient d. which gives the friction 
coefficient for the compressible case with uniform mass flux by 
the relationship 


The detailed steps of the transformations and the derivation 
of Eq. (10) are found in Appendix A of reference 4. 
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A Perturbation Analysis of Shock Flow 
in a Nozzle* 


Roy Gundersen 
Assistant Professor of Mathematics, 
Illinois Institute of Technology, Chicago, Ill. 


May 23, 1959 


B* THE USE of a linearization based on small area variations, 
the motion of a plane shock generated by a uniform com- 
pressive piston motion has been solved for the case of a slowly 
converging or diverging one-dimensional channel.! the 
present paper, these results are extended to the case where the 
channel presents a throat. With E(x) the cross section of the 
tube, it is convenient to let E(x) = Ey + k(x — a), where Ep is 
the original uniform cross section, a is a constant, and k is a small 
parameter. 

At t = 0, a uniform shock is introduced at x = 0, and it is as- 
sumed that the piston is pushed continuously with uniform 
velocity, and the original shock velocity is W. Utilization of the 
general solution for the nonisentropic perturbation of a constant 
state gives the shock-velocity perturbation IV(t). 

= [uecoMikat/(ki + Mike — — 

+ M, + 2m)t?)]/[2ki(1 + Mi?) + 

4Miko — (1) 

where the rest state in front of the shock is denoted by the 

subscript 0 and the perturbed flow behind the shock by the sub- 

script 2, u and ¢ are the fluid and sound velocities, W. = W — 

uo, My) = W/co, = We/co, m = U2/c2, = Wyy(y — 1) with 

c, the specific heat at constant volume, y the adiabatic index 


FORUM 763 
TABLE 1 
1 2 0 0 
4 0.96390 1.0315 1.0702 
8 0.62578 2.0870 3.3350 
10 0.54758 2.5745 4.7016 
TABLE 2 
1 2a 
4 1.8216 a + 0.14973 ka*®/EFy 
8 1.6557 a + 0.13621 ka?/Ey 
1.4359 a + 0.089090 


and k, ko, and k; depend upon the Mach Numbers defined 
above. 
For y = 7/5 and the shock strength defined as the pressure 
ratio P2/P») = 2, Eq. (1) becomes 
W(t) = {(6(2 + — 1)skateo?|/(272? + 652 + 6)} — 
{ [3(s — 1)(2 + 6)(392? — 8z + 18)cotkt?] + 
[14(92? + 31z + 9)(42s + 7)!/2}} 
= Ryco*kat — Rocy*kt? (2) 


which is a quadratic function of ¢ with zeros at ¢ = 0 and ¢,; where 
Cohi/a = 282(92? + 3lz + 9)(42e + 7)!/2/(27s* + + 6) 
(3922 — 8 +18). I(t) is zero at t = 0, increases to a maximum 
at t,/2, and decreases to zero at t;. For larger t, 1 < 0, increas- 
ing in absolute value. Table 1 gives values of 4), Ri, and Re for 
several shock strengths. 

The position of the point where IV’(t) = 0 is given in Table 2. 

For a shock strength less than about 2.2, 17 > 0 throughout 
the entire nozzle. For z > 2.2, the maximum point lies upstream 
of the throat and moves farther upstream for increasing z with a 
limiting position of about 0.72a for a very strong shock. If the 
tube is of constant cross section after the nozzle, the shock will 
asymptotically become uniform. 

For the case a = (—i.e., cross-sectional perturbation E,; = 
kx*—the pressure perturbation directly behind the shock is given 
by 


where 


Ay — 1 + + My + 2m)/(y + 1)?Ke = Mi(1 + + Mo?) + + — 
Mi(y — 1 + 2Mo~*)y~! — + ICL — + — 24+ — y — 14+ (3) 


K, is a monotonic decreasing function of shock strength, and for y = 7/5,0.5 = Ke> 0.361. 
Eq. (3) may be interpreted as a differential relation between shock strength and area, viz., 


dE 
= 1/(2 — 1)Ke(s) = 


— 1) + Rez + — D2 + ¥ + 1) 


Ay — + 1)%2? + — + (y — 14) 


(2 — 1)(Rsz? + Rez + Rs) 
R, = 47? — 6y — 2, 


with R; = 377+ 474+ 1, 


Integration yields Eg(z) = constant where 


{x — 777? — 1097 — 31 
cap 


+ 
4Ri(z + 22y — 


+ 227 — 


For the particular area variation under consideration, this 
relation gives an approximate description of the motion of the 
shock in terms of the cross-sectional area. 


* This work was supported by a grant from the National Science Foun- 
dation. 


val g(z) = (2 — — 1)e + 7 + 1) 


[(y — De + y + + Re + Rs) 
R; = (y+ 1)? 


4(y — 1) 


[R322 + Ry + 
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The Effect of Temperature on Pressures 
Measured in a Hypersonic Wind Tunnel* 


Weston M. Howard 
Aerodynamics Engineer, Boeing Airplane Company, Pilotless 
Aircraft Division, Seattle, Wash. 


June 4, 1959 


— CREEP, a gaseous flow phenomenon which occurs at 
the boundary of a surface and a gas at low pressure, can be of 
considerable importance in hypersonic testing. An important 
effect of thermal creep is the buildup of a pressure gradient 
along the length of a tube. Thus, appreciable error may result 
at high temperatures and low pressures when attempting to 
measure pressures on the surface of a hypersonic-wind-tunnel 
model, since the pressure-lead tubes are generally small and have 
an axial temperature gradient due to high model-surface temper- 
atures. 

Theory for regions of relatively large and small Knudsen 
Numbers (in the measuring tube) has been derived and is well 
established.'~* However, some doubt exists concerning the 
transition region where the Knudsen-Number range is approxi- 
mately 0.1 to 10. For this transition region, a semi-empirical 
equation, based on experiments with hydrogen, was derived by 
Knudsen.‘ > 

The three equations and their regions of applicability are given 
below. 


dp/p = [3.8Kn2/(1 + 4Kn)|\(dT/T), Kn <01 (1) 


dp/p = (1/2){1 + [2.46 (Kn + 3.15)/Kn (Kn + 24.6)]}-? x 
(dT/T), 0.1 < Kn <10 (2) 


dp/p’= (1/2)(dT/T), Kn> 10 (3) 


where K, = Knudsen Number, ratio of mean-free-path length to 
tube inside radius 
p = pressure 
T = temperature 
* The work discussed in this paper was carried out at the Guggenheim 
Aeronautical Laboratory, California Institute of Technology, under the 
sponsorship and with the financial support of the Office, Chief of Ordnance. 
The author wishes to express his appreciation to Dr. H. T. Nagamatsu for 
his guidance throughout this investigation. 
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on the pressure measured through circular tubes. 


In order to check the theory for the range of Knudsen Numbers 
which might be encountered in hypersonic-wind-tunnel testing, 
an experimental program was planned and carried out. The test 
apparatus consisted of stainless-steel tubes of inside radius 1.3 in., 
0.0095 in., and 0.0309 in. The tube lengths were about 11 in. 
Provisions were made for evacuating the tubes and for heating 
one end and cooling the other. The large tube served as a refer- 
ence, since no appreciable pressure gradient would be built up 
for such a large-diameter tube. Therefore, the pressure drop for 
the small-diameter tubes could be found by comparing the small- 
tube cold-end pressure with the large-tube reference pressure. 
A schematic diagram of the apparatus and overall system is shown 
in Fig. 1. 

The results are most clearly described by plotting the data in 
the form 


vs. average Kn, as shown in Fig. 2 where 


( ), refers to hot-end conditions 
( ). refers to cold-end conditions 
( )ave. refers to average conditions 


As can be seen from the plot, the data fall reasonably close to 
the theoretical curves (except for point scatter at large values 
of Kn) even with the rather crude assumption that dp and dT 
can be represented by the large increments p, — p- and T, — T. 
The scatter increases rapidly with Kn due, of course, to the fact 
that the pressure differential is divided by the average pressure, 
which is very small for large Knudsen Numbers. Also, drift and 
outgassing are more detrimental at the low pressures. There- 
fore, correlation with theory is questionable for Kn greater than 
2. However, it appears that the data would follow Knudsen’s 
transition curve up to a somewhat larger Knudsen Number if 
more precise measurements were made. 

Based on these results, it seems reasonable to assume that Eqs. 
(1) and (2) are sufficiently accurate for use in applying temper- 
ature corrections to pressure-test data. The corrections can 
then be made either by using the equations, or by reading the 
ordinate from Fig. 2 (after calculating the Knudsen Number). 
It should be noted that Kn is the average value in the tube; 
therefore, for the case of large corrections, several iterations 
should be made, each time using the new average value of Kn and 
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Acoustical Instability of Jets and Wakes 
for Mach Numbers Above Two 


Erik Mollé-Christensen 
Associate Professor of Aeronautics and Astronautics, 
Massachusetts Institute of Technology, Cambridge, Mass. 


July 8, 1959 


IBNER! AND MILEs? have shown that a discrete interface 

between two moving media may be unstable to supersonic 
disturbances when the Mach Number difference between the 
two media exceeds two. 

This note points out that for a planar ideal jet the acoustical 
instability will lead to sound emission at one frequency f = 
aM/4hV M? — 4 and with wave fronts making an angle with 
the upstream direction of sin—(2/M/). 


INSTABILITY OF A SINGLE INTERFACE 
Consider a planar interface between two streams moving in 
opposite directions with Mach Numbers M, and M2, respec- 
tively, relative to the coordinate system (Fig. 1). A stationary 
disturbance as shown represents a pressure field of the form: 


y > 0; pilx, y) = AF(x + = (1/81) pra. Min(x + a 
y <0; pox, y) = BF(x + Boy) = (1/B2)p2a2Movx(x + Bey) 
where a = velocity of sound, v = vertical velocity, p = density, 
= V/ M2 — 1;7% = 1,2, and A and B are constants. The 
subscripts 1 and 2 refer, respectively, to the upper and lower 
medium. 

The pressure-vertical velocity relation is that given by Ackeret 
(see reference 3) for small disturbances on a steady uniform su- 
personic flow where a linear approximation can be used and no 
viscosity is present. 

The boundary conditions at the interface are: 


pi = praty = 0 (2) 
(zero pressure difference ) 
ve/a,M, = v2/a2Me aty = 0 (3) 


(no difference in streamline slope) 

Substituting Eq. (1) into Eqs. (2) and (3), one obtains: 
-A+B=0 
(2? — B= 0 


The condition for a nontrivial solution for A and B is that the 
determinant of this set of equations vanishes, which yields: 


a2? =1 


Expressing M2 in terms of M4, and M, one can then solve for M, 
in terms of a, d2, pi, p2 and M. For the special case of az = a, and 
p2 = pi, one finds Solution 1: M, = M2: = M/2, and two other 
solutions, of which one or both are real for M > 2. We will 
refer to these as Solutions 2 and 3. It is interesting to note that 
Solution 1 is valid to the second order in streamline slope, since 
the second-order correction is a function of Mach Number and 
streamline slope only, and thus is the same above and below 
the interface, 1, being equal to M2. Solutions 2 and 3 do not 
exist in the second-order approximation. 
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Using the second-order relation between vertical velocity v and 
pressure p, one obtains 
pr = (pai + + — » > 0 
P2 = + + 1)Me* — 482") 
Substituting this into the boundary conditions, one obtains for 
pi = peand ay = az: 

+ [(y + 1) Ms — 482%] /4824) 
The only solution for M@, and M, independent of » is 14, = Mb, 
which gives v1; = %, which corresponds to the previous Solution 2. 
The previous Solutions 2 and 3 must, therefore, be limited to 
much smaller amplitudes than Solution 1, and the latter will 
dominate. The disturbance will thus be a planar wave limited 
in amplitude, due to effects of terms of third and higher orders in 
streamline slope. 

It should also be remarked that these disturbances are only 
unstable for Mach Numbers between M = 2 and M = Das 2, 
as shown by Miles.2, Above M = 2vV/2 they are neutrally stable, 
but the presence of another interface will cause disturbances to 
interact and become unstable also for M > 2V/2. There being 
no length scale in the present arrangement, the wave length and 
wave form of the disturbances are arbitrary. 

Next, introduce a length scale by considering a two-dimensional 


ideal jet between two interfaces at y = Oand y = —A. 


INSTABILITY OF AN IDEAL JET 

On each of the interfaces, the disturbances considered above 
can exist. A disturbance of one layer will now impinge on the 
other layer and create a disturbance, which will be reflected. 
As shown in Ribner’s paper,! this reflection will be infinitely 
strong. Nonlinear effects of third order in streamline slope will, 
however, limit the reflection to some finite amplitude, and a 
periodic disturbance pattern will be formed (Fig. 2). 

The period in x will be \-: 


Az = 4h cot[sin=!(2/M)] = 2hv M? — 4 


With respect to an observer at rest in the medium outside the 
jet, these waves will move in the x-direction with a velocity aM/2, 
and their fundamental frequency w/27 is thus: 

Qwh/2raM = M? — 4 
The angle 6 between the wave-front normal and the downstream 
direction is 

@ = cos~! (2/M) 

These waves will probably not be simply harmonic. The exact 


wave form will have to be determined using nonlinear theory, 
including terms of at least third order in streamline slope. 


APPLICATION TO REAL JETS 
In real jets, the shear layer has finite thickness. The analysis 
given above should be approximately valid when the shear-layer 
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thickness is small compared to the wave length of the disturb- 
ances—i.e., to the jet thickness h. 

With the disturbances there will be associated Reynolds 
stresses in the shear layer; their net effect will be to exert a net 
shearing stress on the layer. The direction of this shearing stress 
can be found by considering the drag on the ripples in the inter- 
face caused by the traveling disturbances. On the outside faces 
of the jet, this drag acts opposite to the direction of flow in the 
jet; it thus will tend to slow down the jet. 
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Generalized Spectral Representation in 
Aeroelasticity 


Arthur K. Cross 
Norair Division, Northrop Corporation, Los Angeles, Calif. ° 
June 4, 1959 


—— a few fragmentary specialized results of the general 
technique treated in this note have appeared in the liter- 
ature,! it is believed that the results of this general method are 
useful enough to justify the brief development and catalog of 
these results presented below. 

This general technique, Generalized Spectral Representation, 
(abbreviated GSR hereafter) is a straightforward generalization 
of the so-called Spectral theorem? (sometimes referred to as the 
Schmidt-Hilbert theorem in integral-equation terminology). It 
will be shown that GSR permits the expansion of matrices asso- 
ciated with various physical systems. The analytical basis and 
applications (with particular reference to aeroelasticity) are 
discussed below. 


ANALYTICAL Basis oF GSR 

Consider the (n X n) matrices G, E, F of rank m where 

G = Et = (E*) = E, F+ =(F*)' = F (1) 
(where the superscripts *, ¢, + denote complex conjugation, 
transposition, and hermitian conjugation, respectively). Then, 
it is well known that m eigenvectors X; exist which satisfy the 
following equations: 

FX; = vjEX, (2) 
= tEX; = 1S Sm,1 (3) 
where, clearly, the hermiticity of E and F guarantee the reality of 


A; and N;. Let the unit rank projection operators be denoted by 
Py 


[;+E/N; = X,Xj+F/\;N;, 


Then, the following properties of the P;’s result from Eqs. (1) to 
(3): 


(4) 


™G)= DY (6) 
i=1 i=1 


(where the ordinary scalar power series 7(A) exists). Hence, 


E = NX FX; 
(7) 


+/X; 
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PRACTICAL APPLICATIONS OF GSR IN AEROELASTICITY 

So-called mass and stiffness matrices are often calculated in 
aeroelastic practice by means of matrices containing data sup- 
plied by tabulated functions or known problem geometry and 
diagonal matrices containing data based on statically determined 
elastic and inertial properties of the system under consideration, 
The latter matrices are pre- and post-multiplied by the former, 
and the resulting mass and stiffness matrices are used to deter- 
mine modal frequencies and deflections. Conversely, Eqs. (7) 
provide a reverse process for determining, for instance, an 
EI or GJ distribution (or ‘‘effective distributions’ for systems 
for which an elastic axis assumption is unrealistic) based on 
vibration data. Such information may be used to determine di- 
rectly the compatibility of observed vibration data and statically 
determined stiffness data. An analogous but less practical pro- 
cedure exists for inertial data. The rest of this note is concerned 
with two brief expositions of the first of these reverse processes 
and a discussion of other regions of application of GSR. 

Although the following results hold quite generally, to fix 
ideas consider a case where the symbols C**, C*’, C’*, C”? denote 
the (m X m) vertical, vertical-roll, roll-vertical, roll, and twist 
submatrices of a (383m X 3m) flexibility matrix C. Numerical 
integration of well-known formulas fcr the elements of these 
matrices provides a sequence of matrix operations by which the 
above matrix sequence may be obtained: W;U(EI)—'W,U'W,, 
W,U(EI)W.V'Ws, 
the ij” of U and V are, respec- 


tively, (y; — 6x; and is the i” wing station 
=1 


numbering from oe zero station in the direction of positive y, 
and the W; are weighting matrices). Furthermore, let h, h’, 
and a denote the (m X 1) (bending, bending slope, and torsion) 
submatrices of the (8m X 1) matrix H. Then, if go is the (3 X 1) 
matrix of zero station deflections, the following matrix formula- 
tion of the equations of motion (in vacuo) defines the (3m X 3) 


matrix Lo: 
H = H — Logo = wtCMH (M = mass matrix) (8) 


The system dynamic loads balance equations define the (3 X 3) 
and (3 X L) matrices L; and Le: 


LiqQo = —1.T 
m m m m t (9) 
T= m,hi, Saini, I,,h;’, > ae 
1 1 1 1 ) 
(T, of course, is an (Z X 1) matrix). 
Finally, if the (ZL X 3m) matrix L; is constructed such that 
L;H = T (10) 


and the matrices 7 and qo are eliminated between Eqs. (8)-(10), 
there results 


= KCMH, K =I — Lo Li + (11) 
Then, according to Eqs. (2) and (7), 
3m 
C= > Ht /H (12) 
i=1 


The diagonal (m X m) matrices (EJ) and (GJ) are then found by 
application of inverses of the known matrices U, Ut, V, V‘t, W; 
to appropriate submatrices of C. It is assumed above that the 
test modes H; have been orthogonalized by one of the standard 
methods—e.g., the Kennedy-Pancu technique—and decoupled 
from the mount frequencies. If measured cantilever modes are 
used, K =J. The required orthogonality of the H; and the H; 
[see Eq. (12)] follows from Eq. (8). 

A completely analogous procedure utilizing uniform beam 
modes is outlined below. First, let v and v be (3 XK m/2) and 


(m/2 X m/2) matrices of which the ij elements are, respectively, 
the i” 
jand the i” 


rigid free vertical deflection shape evaluated at station 
elastic free bending mode shape evaluated at station 
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j. Then, let v' = v', O] and nt = [n‘o, 0, where v and n 
are |((m + 3) X m/2] matrices and mo and n are (3 XK m/2) and 
(m/2 X m/2) torsional counterparts of v and 7. 

Then, the generalized mass matrix B may be calculated as 


follows: 


B = vMvt + vSant + "Sav! — vS,v't — v'S,vt + 
+ — — (13) 


Then let P be the (m X m/2) matrix of which the ij element 
is the i” free test mode shape at station j (the bending com- 
ponent thereof) with an analogous definition holding for the 
torsional component matrix Q which has the same dimensions. 
Furthermore, let Z‘, z‘, and Q denote m XK (m + 3),(m X m) and 
(m X m) matrices, respectively, where 


= (Bee = Bee = BeoBooBoe) { 


and Age = {v"(El)v"t + (15) 


The ee, oe, eo, a0 subscripts represent (m X m),(3 X m),(m X 3), 
(3 X 3) submatrix indices, and, as usual, the primes denote span 
differentiation. 


Since = PMvt + PSan't + QSav' — P’S,v' + ete. (16) 


Bcan be eliminated between Eqs. (13) and (16), and the resulting 
columns 2; of z used to determine A,, in accord with Eqs. (7). 
Then, the diagonal matrices (EJ) and (GJ) are calculated by 
applying inverses of the known v”, v”, n’, n't to appropriate sub- 
matrices of A ¢. 
The usefulness of Eq. (6) for sweeping procedures is obvious if 
k 


it is noted that successive iterates of G — > \:P; converge to 
1 


Me+iXe+1. If G is the unit matrix, which has unity eigenvalues 
and arbitrary eigenvectors, F = E, and Eq. (7) provides a con- 
venient expansion of E~' in terms of any easily constructed set 
of XY; orthogonal with respect to E. In actual calculations, Eqs. 
(6) and (7) are most conveniently utilized after conversion into 
matrix product form. 

It is interesting to note that although there are points of simi- 
larity and, in some special cases, coincidence of the operators 
P; and the “normalized annihilator products’’ of Sylvester,* in 
general these operators differ. Indeed, the GSR of the single 
annihilator product (rank = o if G # /) corresponding to a spec- 
tral multiplicity of order ¢ is given by Eq. (6). In such a case, 
there exist infinitely many sets of ¢ P;’s corresponding to the fact 
that any of the infinity of sets of orthogonal Cartesian bases 
spanning the o dimensional subspace (of the m dimensional con- 
figuration space associated with the metric tensor EZ) associated 
with the above spectral multiplicity forms principal axes of the 
quadric determined by F. 

Simple arguments can be constructed from the latter circum- 
stance which show that whereas close eigenvalues of G are rela- 
tively insensitive to small changes in G, its eigenvectors are 
sensitive in a manner that increases without bound as the roots 
coalesce. The failure to observe (during vibration test) theoretical 
mode shapes associated with two close frequencies should thus not 
be unexpected. The techniques described above can be used to 
determine the important component of error by comparing the 
modal contributions to the stiffness distributions based on vibra- 
tion data with the modal contributions to the assumed distribu- 
tions (calculated using theoretical modes). 
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A Two-Cell Vortex Solution of the 
Navier-Stokes Equations* 


Roger D. Sullivan 
Aeronautical Research Associates of Princeton, Princeton, N.J. 


June 26, 1959 


Breecess" 2 has discussed an exact solution of the Navier- 
Stokes equations for a three-dimensional vortex. In terms 
of a cylindrical coordinate system (7, 6, z) with the radial, tangen- 
tial, and axial velocity components designated by u, v, and w, 
respectively, it can be written 
u=—ar, w=2az, v =(P/2er)[1 (1) 


where a and I are constants (1 being the circulation at infinity), 
and » is the kinematic viscosity. Recently, Rott* * has con- 
sidered the same solution in more detail, including an extension 
to the case of v varying with time, and determination of pressure 
and temperature distributions. 

An analog to Burgers’ solution with two cells has been found in 
the course of a study of a family of solutions in which uw and v are 
functions of r only and w is a function of r multiplied by z. By 
“two cells’’ is meant that the flow does not simply spiral in 
toward the axis and out along it, as shown in Fig. l(a), but has 
a region of reverse flow near the axis, as shown in Fig. 1(b). This 
solution is given by 


“u= —ar+ 6(v/r) [1 e~(ar?/2v) | (2 
w = 2az [1 — 3e~*/2")] (3) 
v = /2nr)[H(ar?/2v)/H(@ )] (4) 


where H is defined by 
x \ t 
H(x) = f, exp +3 (1 —e dt 
The pressure for the two-cell case is 
= po — (p/2){4a%s? + a*r? + 36(v2/r?) X 
r 
[1 (ar? 2v))2} + ef, (v?/r)dr (5) 
where v is given by Eq. (4), while for the one-cell case it is 
b = po — (p/2)[4a2s? + a*r?] + ef, (v?/r)dr 


with v given by Eq. (1). It is to be noted that for both cases the 
axial pressure gradient is 
Op/Os = —4pa*s 
* This research was supported by the USAF through the AFOSR of the 


ARDC, under Contract No. AF 49(638)-255. Reproduction in whole or in 
part is permitted for any purpose of the United States Government. 


Fic. 1. Vortex sketches. (a) Left: one cell. (b) Right: 
two cells. 
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Fic. 2. Velocity distributions. (a) Top: axial velocity. 
(b) Bottom: tangential velocity. 


and that Eqs. (2), (3), (4), and (5) all reduce, for large 7, to the 
corresponding expressions for the one-cell case. 

The distribution of w with ry for the two cases is shown in 
Fig. 2(a), and that of v in Fig. 2(b). A marked reduction in the 
maximum value of v for the two-cell as compared with the one- 
cell case is evident, and there is a large region near the axis where 
v can be considered negligible in the two-cell case in comparison 
with the one-cell case. 

Burgers calculated the dissipation, using only the terms in v, 
and found the result to be of the third order in velocities and 
independent of viscosity. The same is found to be true in the 
two-cell case. 
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On the Accuracy of an Approximate Theory 
for the Ablation of Glassy Materials—A Reply 


Hans A. Bethe and Mac C. Adams 
Avco-Everett Research Laboratory, Everett, Mass. 
June 24, 1959 


T° A RECENT NOTE, Sutton! compared several approximate 
theories with his numerical results for the case of a glass which 
melts and flows without evaporation. He concluded: (1) our 
theory? predicts an ablation rate which is too large by as much 
as 25 per cent, and (2) this same theory will yield significant 
errors when calculating the amount of melt which vaporizes, 
These same conclusions were also quoted in a recent paper by 
Lees.* Unfortunately, Sutton’s comparison was made incor- 
rectly, and both conclusions are wrong. 

To make the comparison, he chose a surface temperature of 
4,000°F. and then calculated the heat-transfer and ablation rates 
for this given temperature. His results indicated that our 
theory led to ablation rates which were too large, by about 25 per 
cent. However, it should be noted the corresponding heat- 
transfer rates were also larger by about 35 per cent. He, there- 
fore, compared ablation rates for two different heat-transfer rates, 
and the comparison is not meaningful. 

The comparison should have been made for equal heat-transfer 
rates, since the heat-transfer rate and not the surface temper- 
ature is the boundary condition for a given flight situation. For 
equal heat-transfer rates, our approximate theory gives an 
ablation rate which is 2 per cent larger than given by Sutton’s 
numerical results. The surface temperature is then 3,660°F. in- 
stead of 4,000°F. This agreement is even closer than should be 
expected, since the appendix? gives a second-order correction of 
4 per cent for this particular glass. Clearly, we are discussing 
small numbers, and Sutton’s numerical results check our theory 
quite adequately. 

Now, regarding the second point, our theory becomes even more 
accurate for situations where some of the melt vaporizes. This 
is true for two reasons. First, the vapor pressure is a strong func- 
tion of temperature, and there is a condition of compatibility 
between the vapor pressure and the gas boundary layer. Asa 
consequence of these conditions, the surface temperature is 
fixed by the vapor pressure and not vice versa, as stated by Sut- 
ton. Second, the error, due to the approximation that the melt 
absorbs energy appropriate to the surface temperature (yn = 1), 
becomes insignificant; with evaporation, this energy constitutes 
only a small fraction of the total energy absorbed. 

We wish to call attention to a report* which compares our 
theory with quartz ablation experiments performed in an air are 
wind tunnel. The experimental data scatter above and below the 
theory by less than 8 per cent. 
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